{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "###### Content under Creative Commons Attribution license CC-BY 4.0, code under BSD 3-Clause License © 2018 parts of this notebook are from [this Jupyter notebook](https://nbviewer.jupyter.org/github/krischer/seismo_live/blob/master/notebooks/Computational%20Seismology/The%20Finite-Difference%20Method/fd_first_derivative.ipynb) by Kristina Garina and Heiner Igel [@heinerigel](https://github.com/heinerigel) which is a supplemenatry material to the book [Computational Seismology: A Practical Introduction](http://www.computational-seismology.org/),  additional modifications by D. Koehn, notebook style sheet by L.A. Barba, N.C. Clementi"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<link href=\"https://fonts.googleapis.com/css?family=Merriweather:300,300i,400,400i,700,700i,900,900i\" rel='stylesheet' >\n",
       "<link href=\"https://fonts.googleapis.com/css?family=Source+Sans+Pro:300,300i,400,400i,700,700i\" rel='stylesheet' >\n",
       "<link href='http://fonts.googleapis.com/css?family=Source+Code+Pro:300,400' rel='stylesheet' >\n",
       "<style>\n",
       "\n",
       "@font-face {\n",
       "    font-family: \"Computer Modern\";\n",
       "    src: url('http://mirrors.ctan.org/fonts/cm-unicode/fonts/otf/cmunss.otf');\n",
       "}\n",
       "\n",
       "\n",
       "#notebook_panel { /* main background */\n",
       "    background: rgb(245,245,245);\n",
       "}\n",
       "\n",
       "div.cell { /* set cell width */\n",
       "    width: 800px;\n",
       "}\n",
       "\n",
       "div #notebook { /* centre the content */\n",
       "    background: #fff; /* white background for content */\n",
       "    width: 1000px;\n",
       "    margin: auto;\n",
       "    padding-left: 0em;\n",
       "}\n",
       "\n",
       "#notebook li { /* More space between bullet points */\n",
       "margin-top:0.5em;\n",
       "}\n",
       "\n",
       "/* draw border around running cells */\n",
       "div.cell.border-box-sizing.code_cell.running { \n",
       "    border: 1px solid #111;\n",
       "}\n",
       "\n",
       "/* Put a solid color box around each cell and its output, visually linking them*/\n",
       "div.cell.code_cell {\n",
       "    background-color: rgb(256,256,256); \n",
       "    border-radius: 0px; \n",
       "    padding: 0.5em;\n",
       "    margin-left:1em;\n",
       "    margin-top: 1em;\n",
       "}\n",
       "\n",
       "\n",
       "div.text_cell_render{\n",
       "    font-family: 'Source Sans Pro', sans-serif;\n",
       "    line-height: 140%;\n",
       "    font-size: 110%;\n",
       "    width:680px;\n",
       "    margin-left:auto;\n",
       "    margin-right:auto;\n",
       "}\n",
       "\n",
       "/* Formatting for header cells */\n",
       ".text_cell_render h1 {\n",
       "    font-family: 'Merriweather', serif;\n",
       "    font-style:regular;\n",
       "    font-weight: bold;    \n",
       "    font-size: 250%;\n",
       "    line-height: 100%;\n",
       "    color: #004065;\n",
       "    margin-bottom: 1em;\n",
       "    margin-top: 0.5em;\n",
       "    display: block;\n",
       "}\t\n",
       ".text_cell_render h2 {\n",
       "    font-family: 'Merriweather', serif;\n",
       "    font-weight: bold; \n",
       "    font-size: 180%;\n",
       "    line-height: 100%;\n",
       "    color: #0096d6;\n",
       "    margin-bottom: 0.5em;\n",
       "    margin-top: 0.5em;\n",
       "    display: block;\n",
       "}\t\n",
       "\n",
       ".text_cell_render h3 {\n",
       "    font-family: 'Merriweather', serif;\n",
       "\tfont-size: 150%;\n",
       "    margin-top:12px;\n",
       "    margin-bottom: 3px;\n",
       "    font-style: regular;\n",
       "    color: #008367;\n",
       "}\n",
       "\n",
       ".text_cell_render h4 {    /*Use this for captions*/\n",
       "    font-family: 'Merriweather', serif;\n",
       "    font-weight: 300; \n",
       "    font-size: 100%;\n",
       "    line-height: 120%;\n",
       "    text-align: left;\n",
       "    width:500px;\n",
       "    margin-top: 1em;\n",
       "    margin-bottom: 2em;\n",
       "    margin-left: 80pt;\n",
       "    font-style: regular;\n",
       "}\n",
       "\n",
       ".text_cell_render h5 {  /*Use this for small titles*/\n",
       "    font-family: 'Source Sans Pro', sans-serif;\n",
       "    font-weight: regular;\n",
       "    font-size: 130%;\n",
       "    color: #e31937;\n",
       "    font-style: italic;\n",
       "    margin-bottom: .5em;\n",
       "    margin-top: 1em;\n",
       "    display: block;\n",
       "}\n",
       "\n",
       ".text_cell_render h6 { /*use this for copyright note*/\n",
       "    font-family: 'Source Code Pro', sans-serif;\n",
       "    font-weight: 300;\n",
       "    font-size: 9pt;\n",
       "    line-height: 100%;\n",
       "    color: grey;\n",
       "    margin-bottom: 1px;\n",
       "    margin-top: 1px;\n",
       "}\n",
       "\n",
       "    .CodeMirror{\n",
       "            font-family: \"Source Code Pro\";\n",
       "\t\t\tfont-size: 90%;\n",
       "    }\n",
       "/*    .prompt{\n",
       "        display: None;\n",
       "    }*/\n",
       "\t\n",
       "    \n",
       "    .warning{\n",
       "        color: rgb( 240, 20, 20 )\n",
       "        }  \n",
       "</style>\n",
       "<script>\n",
       "    MathJax.Hub.Config({\n",
       "                        TeX: {\n",
       "                           extensions: [\"AMSmath.js\"], \n",
       "                           equationNumbers: { autoNumber: \"AMS\", useLabelIds: true}\n",
       "                           },\n",
       "                tex2jax: {\n",
       "                    inlineMath: [ ['$','$'], [\"\\\\(\",\"\\\\)\"] ],\n",
       "                    displayMath: [ ['$$','$$'], [\"\\\\[\",\"\\\\]\"] ]\n",
       "                },\n",
       "                displayAlign: 'center', // Change this to 'center' to center equations.\n",
       "                \"HTML-CSS\": {\n",
       "                    styles: {'.MathJax_Display': {\"margin\": 4}}\n",
       "                }\n",
       "        });\n",
       "</script>\n"
      ],
      "text/plain": [
       "<IPython.core.display.HTML object>"
      ]
     },
     "execution_count": 1,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Execute this cell to load the notebook's style sheet, then ignore it\n",
    "from IPython.core.display import HTML\n",
    "css_file = '../style/custom.css'\n",
    "HTML(open(css_file, \"r\").read())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Finite-difference approximations\n",
    "\n",
    "In the last lecture we covered different approaches to discretize continous media.\n",
    "Now, we only have to understand how to approximate (partial) derivatives by finite-differences."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Approximation of first derivatives\n",
    "\n",
    "The approximation of first derivatives using finite-differences is quite straightforward. Imagine that we have some time-dependent function $f(t)$, for which we want to calculate the first derivatives $\\frac{\\partial f(t)}{\\partial t}$:\n",
    "\n",
    "<img src=\"images/gauss_disc_final.png\" width=\"95%\">\n",
    "\n",
    "First, we discretize $f(t)$ at discrete times\n",
    "\n",
    "\\begin{equation}\n",
    "t = i * dt, \\nonumber\n",
    "\\end{equation}\n",
    "\n",
    "where $dt$ denotes a constant temporal sample interval. By definition, the first derivative of $f(t)$ with respect to $t$ is:\n",
    "\n",
    "\\begin{equation}\n",
    "\\frac{\\partial f(t)}{\\partial t} = \\lim_{dt \\rightarrow 0} \\frac{f(t+dt)-f(t)}{dt}, \\nonumber\n",
    "\\end{equation}\n",
    "\n",
    "In the **finite-difference (FD) approximation**, we neglect the limit $dt \\rightarrow 0$:\n",
    "\n",
    "\\begin{equation}\n",
    "\\frac{\\partial f(t)}{\\partial t} \\approx \\frac{f(t+dt)-f(t)}{dt} \\nonumber\n",
    "\\end{equation}\n",
    "\n",
    "Because $dt$ remains finite, this approximation is called finite-differences. Furthermore, we can distinguish between different finite difference operators, depending on the points involved in the approximation. The above example is a **forward FD operator**\n",
    "\n",
    "\\begin{equation}\n",
    "\\biggl(\\frac{\\partial f(t)}{\\partial t}\\biggr)^+ \\approx \\frac{f(t+dt)-f(t)}{dt}. \\nonumber\n",
    "\\end{equation}\n",
    "\n",
    "Alternatively, we can also define a **backward FD operator**\n",
    "\n",
    "\\begin{equation}\n",
    "\\biggl(\\frac{\\partial f(t)}{\\partial t}\\biggr)^- \\approx \\frac{f(t)-f(t-dt)}{dt}. \\nonumber\n",
    "\\end{equation}\n",
    "\n",
    "By taking the arithmetic average of the forward and backward operator\n",
    "\n",
    "\\begin{equation}\n",
    "\\biggl(\\frac{\\partial f(t)}{\\partial t}\\biggr) \\approx \\frac{\\left(\\frac{\\partial f(t)}{\\partial t}\\right)^- + \\left(\\frac{\\partial f(t)}{\\partial t}\\right)^+}{2},\\nonumber\n",
    "\\end{equation}\n",
    "\n",
    "we get the **central FD operator**\n",
    "\n",
    "\\begin{equation}\n",
    "\\biggl(\\frac{\\partial f(t)}{\\partial t}\\biggr) \\approx \\frac{f(t+dt)-f(t-dt)}{2dt}. \\nonumber\n",
    "\\end{equation}\n",
    "\n",
    "##### Exercise\n",
    "\n",
    "Prove for the quadratic function\n",
    "\n",
    "\\begin{equation} \n",
    "g(t) = b t^2 \\nonumber\n",
    "\\end{equation} \n",
    "\n",
    "where $b$ is a constant time-independent parameter, that in the limit $dt \\rightarrow 0$, the forward, backward and central FD operators lead to the correct temporal first derivative:\n",
    "\n",
    "\\begin{equation} \n",
    "\\frac{\\partial g(t)}{\\partial t} = 2 b t \\nonumber\n",
    "\\end{equation} "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "code_folding": [
     0
    ]
   },
   "outputs": [],
   "source": [
    "# Import Libraries\n",
    "import numpy as np\n",
    "import math\n",
    "import matplotlib.pyplot as plt\n",
    "from pylab import rcParams"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As an example, let's try to calculate the first temporal derivative of the Gaussian function\n",
    "\n",
    "\\begin{equation} \n",
    "f(t)=\\dfrac{1}{\\sqrt{2 \\pi a}}e^{-\\dfrac{(t-t_0)^2}{2a}},\\nonumber\n",
    "\\end{equation} \n",
    "\n",
    "where $a$ denotes the half-width of the Gaussian and $t_0$ a time-shift of the maximum."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "code_folding": []
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAt0AAAFNCAYAAADcudMsAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJzs3Xd8VFX+//HXJx0SOqEFSKFHUJogIBCs2BXL2tFVsK676za3qLu6rq7uqmsHFbsi6qrooohgKFKkCEgnJJTQOySE1PP7I8P3F2OAALm5k8z7+XjMg7n3njvzHrxOPpyce4455xAREREREe+E+R1ARERERKS2U9EtIiIiIuIxFd0iIiIiIh5T0S0iIiIi4jEV3SIiIiIiHlPRLSIiIiLiMRXdIiI1mJl9YWbDPXrtv5vZDjPb4sXrH+F9XzKz+6vzPUVEvGaap1tEpPLM7Grg10BXIBfIAt4AXnS16AvVzNoAq4BE59w2D9/nJuBW59zpXr2HiEgwUE+3iEglmdlvgP8ATwAtgObA7cAAIMrHaF5IBHZ6WXCLiIQSFd0iIpVgZg2Ah4A7nXMfOuf2u1LfO+euc87lB9pdYGbfm9k+M9tgZn8t8xppZpZd7nXXmtlZged9zGxe4NytZvZkYH+Mmb1tZjvNbI+ZzTWz5oFj6WZ2a+B5OzObEmi3w8zeMbOG5d7rt2a22Mz2mtn7ZhZTwWc9C5gEtDKzHDN7vRLZ/2pm48zsTTPbb2ZLzax3mbZtzOy/ZrY9kO85M+sCvAT0C7zPnkDb183s72XOHWFmGWa2y8zGm1mrMsecmd1uZqvNbLeZPW9mdkz/cUVEqoGKbhGRyukHRAOfHqVdLnAj0BC4ALjDzC6t5Hv8B/iPc64+0A4YF9g/HGgAtAGaUNq7nlfB+QY8CrQCugTa/7Vcm6uAoUAycDJwU/kXcc59DZwHbHLOxTnnftLmMC4GxlL62ccDzwGYWTjwObAOSAISgLHOueWBzzIr8D4Ny7+gmZ0R+ExXAS0DrzG2XLMLgVOBUwLtzq1kXhGRaqOiW0SkcpoCO5xzRYd2mNnMQM9znpkNAnDOpTvnfnDOlTjnFgPvAYMr+R6FQHsza+qcy3HOzS6zvwnQ3jlX7Jyb75zbV/5k51yGc26Scy7fObcdeLKC937GObfJObcL+Azofgx/B0czwzk3wTlXDLxFaREM0IfSfwj8zjmX65w76JybUcnXvA4Y45xbEPhtwh8p7RlPKtPmMefcHufceuAbqvYziYhUCRXdIiKVsxNoamYRh3Y45/oHemd3Evg+NbO+ZvZNYBjFXkp7cptW8j1uAToCKwJDSC4M7H8LmAiMNbNNZva4mUWWP9nMmpnZWDPbaGb7gLcreO+yM5EcAOIqma0yyr92TODvqw2wruw/WI5BK0p7twFwzuVQ+vedcIT3rcrPJCJSJVR0i4hUziwgH7jkKO3epXRoRRvnXANKxywfGmOcC9Q91DAw7CL+0LZzbrVz7hqgGfBP4EMzi3XOFTrn/uacSwX6Uzqc4sYK3vtRwAEnB4aoXF/mvU/UEbMfxQagbdl/sJRxtBlfNlF6U+eh942ltNd/YyXfW0QkKKjoFhGpBOfcHuBvwAtmdoWZxZlZmJl1B2LLNK0H7HLOHTSzPsC1ZY6torT394JAT/VfKB0nDoCZXW9m8c65EmBPYHexmQ0xs26BQncfpcNNiiuIWQ/IAfaYWQLwuyr58JXIfhTfAZuBx8wsNnBj6IDAsa1AazM73Owv7wI3m1l3M4sG/gHMcc6tPe5PIiLiAxXdIiKV5Jx7HLgX+D2wjdKCcRTwB2BmoNmdwENmth94gP9/MyTOub2B469Q2lObC5SdEWQosNTMcii9qfJq59xBSqcn/JDSgns5MJXSoSPl/Q3oCewF/gf894Q/dOWzH+ncYuAioD2wPnDezwKHpwBLgS1mtqOCcycD9wMfUVq4twOuPpHPIiLiBy2OIyIiIiLiMfV0i4iIiIh4TEW3iIiIiIjHVHSLiIiIiHhMRbeIiIiIiMdUdIuIiIiIeKyihQpqpIYNG7r27dv7HUOCTG5uLrGxsUdvKCFF14VURNeFVETXhVRk/vz5O5xzlV0gDKhFRXfz5s2ZN2+e3zEkyKSnp5OWluZ3DAkyui6kIroupCK6LqQiZrbuWM/R8BIREREREY+p6BYRERER8ZiKbhERERERj6noFhERERHxmIpuERERERGPqegWEREREfGYp0W3mQ01s5VmlmFm9x2h3RVm5sysd5l9fwyct9LMzvUyp4iIiIiIlzybp9vMwoHngbOBbGCumY13zi0r164ecA8wp8y+VOBq4CSgFfC1mXV0zhV7lVdERERExCte9nT3ATKcc5nOuQJgLHBJBe0eBh4HDpbZdwkw1jmX75zLAjICryciIiIiUuN4uSJlArChzHY20LdsAzPrAbRxzn1uZr8td+7scucmeBVURCRU5BcVs2ZbLpv35rH7QCElzlEnMpz4etGkNI0lvl40ZuZ3TBGRWsfLoruib233fwfNwoCngJuO9dwyrzESGAkQHx9Penr68eSUWiwnJ0fXhfxEqF0X6/cVM39rMYt3FLN+XwnFP/k2/f/qRxmpTcI4OT6CXs3CiY4InQI81K4LqRxdF1JVvCy6s4E2ZbZbA5vKbNcDugLpgV6VFsB4M7u4EucC4JwbDYwG6NSpk0tLS6vC+FIbpKeno+tCyguF66KouIRPF27izVlrWZSdixn0bNuIoT0ak9qyPm0a16VR3UjCzDhYWMyWfQfJ3J7LgvW7mblmJ7M35xMbFc6lPRK4bVA72jap6/dH8lwoXBdy7HRdSFXxsuieC3Qws2RgI6U3Rl576KBzbi/Q9NC2maUDv3XOzTOzPOBdM3uS0hspOwDfeZhVRKRWcM7x2eLNPDVpFVk7cunQLI4HL0rlku4JNI6NOux5HZrXY2CHeIb3T6KkxPHd2l18OD+bD+ZlM3buBob1SOB3QzvRrF5MNX4aEZHaw7Oi2zlXZGZ3AxOBcGCMc26pmT0EzHPOjT/CuUvNbBywDCgC7tLMJSIiR5a5PYf7P13Ctxk76dyiHqNv6MXZqc2PeYx2WJhxWkoTTktpwu/O7cSoqZm8NXstXyzZwr1nd+Sm/kmEhYXOsBMRkargZU83zrkJwIRy+x44TNu0ctuPAI94Fk5EpJZwzvHBvGweGL+EyPAwHr60K9f2aUt4FRTGzevH8MBFqVx/Wlse+nwZD32+jK+WbeFfV55C60a1f8iJiEhV0YqUIiI1WE5+Eb9+fyG//2gxvRIbMfnewdxwWmKVFNxlpcTH8dpNp/L45SfzQ/ZeLnx2Bt9m7KjS9xARqc1UdIuI1FCb9+ZxxYszGb9oE785uyNv/rwvzep7N+bazLjq1Db8756BxMdFc+OY7xgzIwvnjjAdioiIACq6RURqpOWb93HZ8zPJ3p3H6zf34Rdndqjy3u3DSWoay8d3DeCMzs146PNlPDh+KSUlKrxFRI5ERbeISA0zb+0urnppFg7HuNv6MahjfLVniIuOYNT1vbj19GTenLWO3364iKLikmrPISJSU3h6I6WIiFSt+et2MXzMdzSvH8Pbt/alVcM6vmUJCzP+fEEX6teJ5MlJqziQX8yz1/YgMlz9OSIi5embUUSkhliwfjfDx8ylWf0Y3ht5mq8F9yFmxj1nduD+C1P5cukWfjNukYaaiIhUQD3dIiI1wKqt+xk+5juaxkXx3ojTaO7hDZPH45bTk8kvKubxL1fSoE4kD11y0jHPDy4iUpup6BYRCXJb9x3kpjHfUScynLdv7UuLBsFVcB9yx+B27D1QyKhpmTSqG8m953TyO5KISNBQ0S0iEsRy8ou4+bW57M0r5P3b+gX1gjRmxn3ndWZXbgHPTMkgOT6Wy3q09juWiEhQ0JhuEZEgVVLiuOe971m5dT8vXN+LrgkN/I50VGbGI5d1o29yY/7w0Q8sWL/b70giIkFBRbeISJB6evJqpqzYxoMXpTLYh2kBj1dURBgvXd+LFvVjGPnmfDbtyfM7koiI71R0i4gEoUnLtvLM5NVc0as1N5yW6HecY9YoNopXh/fmYGExd76zgIIizeEtIqFNRbeISJBZuyOXe99fSLeEBvz90q41dhaQDs3r8fgVJ7Nwwx4e+2KF33FERHyloltEJIgUFJVwz9jvCQ83XrqhFzGR4X5HOiHnd2vJTf2TGPNtFl8u2ex3HBER36joFhEJIk9OWsXi7L08NuxkEoJg8Zuq8Kfzu3BKm4b87oPFbNh1wO84IiK+UNEtIhIkvs3Ywahpa7i2b1uGdm3hd5wqExURxnPX9MABv/1AK1aKSGhS0S0iEgR25xZw77iFpDSN5f4LUv2OU+XaNK7LgxelMidrF6/OyPI7johItVPRLSISBP708Q/szi3kmWt6UCeqZo/jPpwrerXmnNTmPDFxJSu37Pc7johItVLRLSLiswk/bOaLJVv41dkdOKlV8C+Ac7zMjEeHdaN+nQh+9f5CTSMoIiFFRbeIiI925xbwwKdL6JbQgJEDU/yO47kmcdE8Ouxklm/ex3NTVvsdR0Sk2qjoFhHx0UOfL2NvXiFPXHkyEeGh8ZV8dmpzhvVI4MWpa1i1VcNMRCQ0hMY3vIhIEJq8fCsff7+RO9Pa07lFfb/jVKs/X9CFuOgI/vjfHzSbiYiEBE+LbjMbamYrzSzDzO6r4PjtZvaDmS00sxlmlhrYn2RmeYH9C83sJS9ziohUt5z8Iv788RI6Na/HXUPa+x2n2jWJi+YvF6Qyf91u3v1uvd9xREQ851nRbWbhwPPAeUAqcM2horqMd51z3Zxz3YHHgSfLHFvjnOseeNzuVU4RET88PWkVW/cf5LHLuxEVEZq/dBzWM4EB7Zvwzy9WsGXvQb/jiIh4ystv+j5AhnMu0zlXAIwFLinbwDm3r8xmLKDfMYpIrbdiyz5em7mWq09tQ4+2jfyO4xsz45FLu1FQXMJfxy/1O46IiKe8LLoTgA1ltrMD+37EzO4yszWU9nTfU+ZQspl9b2ZTzWyghzlFRKqNc44HPllK/ZgIfn9uZ7/j+C6paSz3nNmBL5duYfrq7X7HERHxjDnnTeeymV0JnOucuzWwfQPQxzn3i8O0vzbQfriZRQNxzrmdZtYL+AQ4qVzPOGY2EhgJEB8f32vcuHGefBapuXJycoiLi/M7hgQZP6+LbzcW8vIPBdx8UhSD20T6kiHYFJY4/jwjj3CDhwfUISLMfMmh7wupiK4LqciQIUPmO+d6H8s5EV6FobRnu02Z7dbApiO0Hwu8COCcywfyA8/nB3rCOwLzyp7gnBsNjAbo1KmTS0tLq6rsUkukp6ej60LK8+u62JtXyG//nU73Ng25/7r+hPlUXAajsJZbueWNeayNTORWn+Yr1/eFVETXhVQVL4eXzAU6mFmymUUBVwPjyzYwsw5lNi8AVgf2xwduxMTMUoAOQKaHWUVEPPfs5NXszC3g75d2VcFdzhmdm5HWKZ6nv17Ntv26qVJEah/Pim7nXBFwNzARWA6Mc84tNbOHzOziQLO7zWypmS0E7gWGB/YPAhab2SLgQ+B259wur7KKiHht7Y5c3pi1lqt6taFrQu1d6v14mRkPXJhKflExj3+50u84IiJVzsvhJTjnJgATyu17oMzzXx7mvI+Aj7zMJiJSnR77YgWR4WH85pyOfkcJWinxcfz89GRGTc3khtMSOaVNQ78jiYhUmdCcHFZEpBrNydzJl0u3cMfgdjSrH+N3nKD2izM60CQ2in9MWI5XN/qLiPhBRbeIiIdKShx//99yWjaI8e0GwZokLjqCX53VgTlZu5i8fJvfcUREqoyKbhERD32ycCM/bNzL74d2ok5UuN9xaoSr+7QlpWksj36xnKLiEr/jiIhUCRXdIiIeySsovSnw5NYNuOSUn6wNJocRGR7G74d2Zs32XMbNy/Y7johIlVDRLSLikTHfZrFl30H+ckGqpgg8Ruee1JzeiY14ctIqcvOL/I4jInLCVHSLiHhgz4ECXpq6hrO6NKNPcmO/49Q4ZsafLujCjpx8Rk/TMg0iUvOp6BYR8cBLUzPJyS/it+d28jtKjdWzbSPO79aC0dMytWCOiNR4KrpFRKrYtn0HeX1mFpec0orOLer7HadG+/25nSkoLuGFb9b4HUVE5ISo6BYRqWLPTFlNUbHj12drIZwTldQ0lit7tebdOevZuCfP7zgiIsdNRbeISBVatzOXsd9t4GentiGxSazfcWqFX5zZAYBnJ6/2OYmIyPFT0S0iUoWe/no1EeHGPYFCUU5cQsM6XNu3LR/Mz2btjly/44iIHBcV3SIiVWTllv18snAjw/sn0VzLvVepO4e0IzLcePrrVX5HERE5Liq6RUSqyNNfryIuKoI7BrfzO0qt06xeDMP7J/Hpok2s2rrf7zgiIsdMRbeISBVYsWUfXyzZws0DkmhYN8rvOLXS7YPaERsVwZNfqbdbRGoeFd0iIlXgmcmrqRcdwc9PT/Y7Sq3VKDaKW05P5sulW1iyca/fcUREjomKbhGRE7Riyz4m/LCFm9TL7blbBiZTPyaCZ6doJhMRqVlUdIuInKBnJ2cQFx3BLerl9lz9mEhuHpDMxKVbWb55n99xREQqTUW3iMgJWLllP//7YTM39Vcvd3X5+YBk4qIjeO6bDL+jiIhUmopuEZET8MyU1erlrmYN6kYyvH8iE37YzGrNZCIiNYSKbhGR47Ryy34mBHq5G8Wql7s63XJ6CnUiw9XbLSI1hopuEZHj9MyU1dSNDFcvtw8ax0ZxQ79EPlu0icztOX7HERE5KhXdIiLHIWNbaS/3cPVy+2bEwBSiIsJ4/ps1fkcRETkqT4tuMxtqZivNLMPM7qvg+O1m9oOZLTSzGWaWWubYHwPnrTSzc73MKSJyrF5MzyQmIpxbB6b4HSVkNY2L5rq+iXyycCPrdx7wO46IyBF5VnSbWTjwPHAekApcU7aoDnjXOdfNOdcdeBx4MnBuKnA1cBIwFHgh8HoiIr7L3n2ATxdu5Oo+bWisXm5f3TYohfAw44V0je0WkeDmZU93HyDDOZfpnCsAxgKXlG3gnCs7yWos4ALPLwHGOufynXNZQEbg9UREfPfytEzMSoc3iL+a1Y/hmlPb8OH8bDbtyfM7jojIYXlZdCcAG8psZwf2/YiZ3WVmayjt6b7nWM4VEaluO3LyGTt3A5f1SKBVwzp+xxFgxKDSf/y8Mj3L5yQiIocX4eFrWwX73E92OPc88LyZXQv8BRhe2XPNbCQwEiA+Pp709PQTySu1UE5Ojq4L+YkTuS4+XFVAQVEJPevs1LUVRPq2COft2Vn0jN5KXFRFP0KOTt8XUhFdF1JVvCy6s4E2ZbZbA5uO0H4s8OKxnOucGw2MBujUqZNLS0s7gbhSG6Wnp6PrQso73uti38FCfvHNFM7v1pKrL+hZ9cHkuLXqsp9znppGRlgCv0rreFyvoe8LqYiuC6kqXg4vmQt0MLNkM4ui9MbI8WUbmFmHMpsXAKsDz8cDV5tZtJklAx2A7zzMKiJyVG/NWsf+/CLuSGvndxQpp2PzepzVpTmvz1xLbn6R33FERH7Cs6LbOVcE3A1MBJYD45xzS83sITO7ONDsbjNbamYLgXspHVqCc24pMA5YBnwJ3OWcK/Yqq4jI0eQVFDNmRhaDO8bTNaGB33GkAncOaceeA4WMnbvh6I1FRKqZl8NLcM5NACaU2/dAmee/PMK5jwCPeJdORKTyxs3bwM7cAu4a0t7vKHIYPds2om9yY16ZnskNpyUSFaH130QkeOgbSUTkKAqLSxg9LZNTkxrRJ7mx33HkCO5Ia8fmvQf5dOFGv6OIiPyIim4RkaMYv3ATG/fkcWeaermD3eCO8aS2rM9LU9dQUvKTSa9ERHyjoltE5Aicc7w8PZPOLeqR1ine7zhyFGbGHWntWLM9l6+WbfU7jojI/1HRLSJyBNNX72DFlv3cOjAFs+Ob/1mq13ldW5DYpC4vTl2Dc+rtFpHgoKJbROQIXp6eSfP60Vx8Siu/o0glRYSHcdugdizasIdZmTv9jiMiAqjoFhE5rOWb9zF99Q6G90/STBg1zLCeCTSNi2L0tEy/o4iIACq6RUQO6+XpmdSNCue6Pol+R5FjFBMZzvB+SaSv3M7KLfv9jiMioqJbRKQiW/YeZPzCTVzVuw0N6kb6HUeOw/WnJVInMpxXpqu3W0T8p6JbRKQCr89cS4lz3HJ6st9R5Dg1io3iqt6t+WThRrbtO+h3HBEJcSq6RUTKyckv4p056ziva0vaNK7rdxw5AT8/PZniEsdrM9f6HUVEQpyKbhGRcsbN3cD+g0XcOlC93DVdYpNYhnZtwTuz15GTX+R3HBEJYSq6RUTKKCou4dUZWZya1IgebRv5HUeqwIiBKew7WMS4uRv8jiIiIUxFt4hIGV8u3cLGPXmMGJjidxSpIj3aNqJPUmNenZFFUXGJ33FEJESp6BYRCXDO8fK0TJKbxnJWl+Z+x5EqNGJQChv35DFhyRa/o4hIiFLRLSISMHftbhZl7+WW05MJC9OS77XJmZ2bkRIfy+hpWhpeRPyholtEJGD0tEwa1Y3k8p6t/Y4iVSwszBgxMIUlG/cxO3OX33FEJASp6BYRAdZsz+Hr5Vu5oV8SdaLC/Y4jHrisx6Gl4df4HUVEQpCKbhER4NUZWURFhHFjPy35XlvFRIZzY78kvlm5nVVbtTS8iFQvFd0iEvJ25uTz0fxsLu+ZQNO4aL/jiIeuPy2RmMgwLQ0vItVORbeIhLy3Zq8jv6iEW07XNIG1XePYKK7q3YZPvt+kpeFFpFqp6BaRkHawsJg3Z63jzM7NaN8szu84Ug1uOT2ZwpISXtfS8CJSjVR0i0hI+++CjezKLWDEIPVyh4rEJrEMPakF78xZT66WhheRauJp0W1mQ81spZllmNl9FRy/18yWmdliM5tsZolljhWb2cLAY7yXOUUkNJWUOF6Znkm3hAb0TW7sdxypRrcOTGFvXiEfzNPS8CJSPTwrus0sHHgeOA9IBa4xs9Ryzb4HejvnTgY+BB4vcyzPOdc98LjYq5wiErqmrNhG5o5cRgxKwUyL4YSSXomN6JXYiFe/zaK4RIvliIj3vOzp7gNkOOcynXMFwFjgkrINnHPfOOcOBDZnA1qRQkSqzejpmSQ0rMP5XVv4HUV8MGJgMht25TFxqZaGFxHveVl0JwBlf2+XHdh3OLcAX5TZjjGzeWY228wu9SKgiISuRRv28F3WLm4ekEREuG5vCUVnp7YgsUldRk/L1NLwIuK5CA9fu6Lf1Vb4rWZm1wO9gcFldrd1zm0ysxRgipn94JxbU+68kcBIgPj4eNLT06skuNQeOTk5ui7kJ3Jycnjho9nUiYCE/HWkp6/3O5L4ZGCzIt5efoBXPplCy8g8fV/IT+jniFQVL4vubKBNme3WwKbyjczsLODPwGDnXP6h/c65TYE/M80sHegB/Kjods6NBkYDdOrUyaWlpVXtJ5AaLz09HV0XUt4HE6Ywb2seIwamcN5ZXfyOIz7qU1DE549NYV5OA65pE67vC/kJ/RyRquLl71TnAh3MLNnMooCrgR/NQmJmPYBRwMXOuW1l9jcys+jA86bAAGCZh1lFJIRMWldImBk3DUjyO4r4rG5UBNf3TeSrZVvZklvidxwRqcU8K7qdc0XA3cBEYDkwzjm31MweMrNDs5E8AcQBH5SbGrALMM/MFgHfAI8551R0i8gJ25tXyLTsIi46pRUtG9TxO44EgRv7JxIZFsZX6wr9jiIitZiXw0twzk0AJpTb90CZ52cd5ryZQDcvs4lIaHrvu/UcLIZbByb7HUWCRLN6MVzaoxWfLMhmd24BjWKj/I4kIrWQbtkXkZBRUFTCa99mkdokjJNaNfA7jgSRWwemUFACb89e53cUEamlVHSLSMj4fPEmtu7LZ2hSpN9RJMh0bF6Pk5uG88astRwsLPY7jojUQiq6RSQkOOcYPS2TDs3i6NY03O84EoSGJkeyI6eATxdu9DuKiNRCRy26zayumd1vZi8HtjuY2YXeRxMRqTozMnawYst+Lfkuh9WlcRipLevz8vQsSrQ0vIhUscr0dL8G5AP9AtvZwN89SyQi4oGXp2cRXy+aS7q38juKBCkzY8SgZDK25TB11Xa/44hILVOZorudc+5xoBDAOZdHxatNiogEpRVb9jFt1XZu6p9EdISGlsjhXXhyK1rUj+Hl6Zl+RxGRWqYyRXeBmdUhsIS7mbWjtOdbRKRGeGV6FnUiw7mub1u/o0iQiwwP4+YBScxcs5MlG/f6HUdEapHKFN0PAl8CbczsHWAy8HtPU4mIVJGt+w7y6cKNXNW7NQ3rav5lObpr+rYlLjqCV9TbLSJV6KhFt3NuEjAMuAl4D+jtnEv3NpaISNV4Y+ZaikocPz9di+FI5dSPieRnp7bhs8Wb2bQnz+84IlJLHLboNrOehx5AIrAZ2AS0DewTEQlquflFvD17HeemtiCxSazfcaQGuXlAEgCvz1zraw4RqT2OtAz8vwN/xgC9gUWU3kB5MjAHON3baCIiJ2bcvA3sO1jEyMEpfkeRGqZ1o7qc360l781Zzy/OaE+9GC2oJCIn5rA93c65Ic65IcA6oKdzrrdzrhfQA8ioroAiIsejqLiEV2dk0TuxET3bNvI7jtRAIwYmsz+/iPfnbvA7iojUApW5kbKzc+6HQxvOuSVAd+8iiYicuC+WbCF7dx4jBqmXW47Pya0b0je5MWNmZFFYXOJ3HBGp4SpTdC83s1fMLM3MBgdWplzudTARkeN1aMn35KaxnN2lud9xpAYbOSiFTXsPMuGHzX5HEZEarjJF983AUuCXwK+AZYF9IiJBaU7WLn7YuJdbByYTFqa1vOT4DenUjJT4WF6enolzWhpeRI5fZaYMPOice8o5d1ng8ZRz7mB1hBMROR6jp2XSJDaKy3u29juK1HBhYcaIgSks2biP2Zm7/I4jIjXYUYtuM8sys8zyj+oIJyJyrFZv3c+UFdu4sV8SMZFa8l1O3GU9EmgSG6XFckTkhBxpysBDepd5HgNcCTT2Jo6IyIl5ZXoW0RFh3NBMtp5kAAAgAElEQVQv0e8oUkvERIZzQ79Env56NRnb9tO+WT2/I4lIDVSZ4SU7yzw2OueeBs6ohmwiIsdk276DfPz9Rq7s3ZrGsVryXarODaclEh0RxqszsvyOIiI1VGWGl/Qs8+htZrcD+me+iASdN2atpbCkhFtO1zSBUrWaxEVzea/WfLRgI9v35/sdR0RqoMrMXvLvMo9HgZ7AVV6GEhE5VqVLvq/n3NQWJDfVku9S9W45PZnC4hLemr3O7ygiUgNVZkz3Lc65H909YmbJHuURETkuH8zbwN68Qi2GI55pFx/HmZ2b89astdwxuB11onSjrohUXmV6uj+s5L6fMLOhZrbSzDLM7L4Kjt9rZsvMbLGZTTazxDLHhpvZ6sBjeGXeT0RCU1FxCa/MyKJXYiN6JWrJd/HOyEEp7D5QyEcLsv2OIiI1zGF7us2sM3AS0MDMhpU5VJ/SWUyOyMzCgeeBs4FsYK6ZjXfOLSvT7Hugt3PugJndATwO/MzMGgMPUjpzigPmB87dfWwfT0RCwZdLS5d8v//CVL+jSC13alIjTmndgFdnZHFtn7ZafElEKu1IPd2dgAuBhsBFZR49gRGVeO0+QIZzLtM5VwCMBS4p28A5941z7kBgczZwaCWLc4FJzrldgUJ7EjC0ch9JREKJc46XA0u+n6Ul38VjZsaIQSlk7cjl6+Vb/Y4jIjXIYXu6nXOfAp+aWT/n3KzjeO0EYEOZ7Wyg7xHa3wJ8cYRzE44jg4jUcnOydrEoey9/v7Qr4ep1lGow9KQWJDSswyvTszjnpBZ+xxGRGuJIw0t+75x7HLjWzK4pf9w5d89RXruin37uMO91PaVDSQYfy7lmNhIYCRAfH096evpRIkmoycnJ0XVRyz01/yD1IqFZbibp6ZWbQ1nXhVTkWK6LwS2KeXfFLsZ8MpmUhrqhsjbT94VUlSPNXrI88Oe843ztbKBNme3WwKbyjczsLODPwGDnXH6Zc9PKnZte/lzn3GhgNECnTp1cWlpa+SYS4tLT09F1UXut3LKfRV9O41dndeCcMztW+jxdF1KRY7kueucX8dmjk5l/oBE/v7Snt8HEV/q+kKpypOElnwX+fOM4X3su0CEwveBG4Grg2rINzKwHMAoY6pzbVubQROAfZnZoGoJzgD8eZw4RqaVGTV1DnchwhvdL8juKhJi46Aiu65vI6GlrWL/zAG2b1PU7kogEuSMNL/mMwwwHAXDOXXykF3bOFZnZ3ZQW0OHAGOfcUjN7CJjnnBsPPAHEAR+YGcB659zFzrldZvYwpYU7wEPOuV3H8sFEpHbL3n2ATxdtYni/JBppyXfxwc0DkhgzI4tR09bwyGXd/I4jIkHuSMNL/nWiL+6cmwBMKLfvgTLPzzrCuWOAMSeaQURqp5enZRJmMGKQ1uoSfzSvH8PlvVrzwfxsfnlWB5rVO+psuiISwg47ZaBzbuqhBzAL2A3sAmYF9omI+GJHTj5j527g0u4JtGxQx+84EsJuG5RCUXEJY2as9TuKiAS5o65IaWYXAGuAZ4DngAwzO8/rYCIih/PGzLUUFJdw2+B2fkeREJfUNJbzu7Xk7dnr2JtX6HccEQlilVkG/t/AEOdcmnNuMDAEeMrbWCIiFcvJL+KNmWs5N7UF7ZvF+R1HhDvS2pGTX8Tbs9f5HUVEglhliu5tzrmMMtuZwLbDNRYR8dK7c9ax72ARt6epl1uCw0mtGpDWKZ4xM7LIKyj2O46IBKnKFN1LzWyCmd1kZsOBz4C5ZjbMzIZ5nE9E5P/kFxXzyvQs+rdrQvc2Df2OI/J/7kxrz87cAsbN23D0xiISkipTdMcAWyldLTIN2A40Bi4CLvQsmYhIOR8v2Mi2/fncmdbe7ygiP3JqUiN6JTZi9LRMCotL/I4jIkHoSFMGAuCcu7k6goiIHElxiWPUtEy6JTRgQPsmfscR+REz4860dtzyxjw+W7SJYT1b+x1JRILMUYvuwIqSvwCSyrY/2uI4IiJV6cslW8jakcsL1/UksJiWSFA5o3MzOreox4vpa7i0ewJhYbpOReT/q8zwkk+AtcCzlM5kcughIlItnHO8ODWDlKaxnHtSC7/jiFTIzLgjrR2rt+Xw9fKtfscRkSBTmaL7oHPuGefcN+UWzBERqRbpK7ezZOM+bh/cjnD1HkoQu6BbS9o0rsML6WtwzvkdR0SCSGWK7v+Y2YNm1s/Meh56eJ5MRITSXu5npqwmoWEdLuuZ4HcckSOKCA9j5KB2LNywh5lrdvodR0SCyFHHdAPdgBuAM4BDt2S7wLaIiKe+zdjJ9+v38PdLuxIZXpl+AhF/XdmrNc9NWc0zk1czoH1Tv+OISJCoTNF9GZDinCvwOoyISHnPTF5Ni/oxXNlbs0FIzRATGc4dg9vx18+WMTtzJ6elaLYdEanc8JJFgFahEJFqNztzJ9+t3cXtg1OIjgj3O45IpV3dpy3x9aJ5ZvJqv6OISJCoTNHdHFhhZhPNbHzg8anXwUREnp2ymvh60Vzdp63fUUSOSUxkOLcNSmHmmp3MXbvL7zgiEgQqU3Q/SOkQk38ATwLfAVoOTkQ8NX/dLr7N2Mltg1KIiVQvt9Q81/VNpGlclHq7RQSoRNEdmB5wL3AB8DpwJvCSt7FEJNQ9MzmDxrFRXNtXvdxSM9WJCmfEwBSmr97BgvW7/Y4jIj47bNFtZh3N7AEzWw48B2wAzDk3xDn3bLUlFJGQs2jDHqau2s6IgSnUjarM/d4iwen60xJpVDeSZ9XbLRLyjtTTvYLSXu2LnHOnBwrt4uqJJSKh7Nkpq2lYN5Ib+iX6HUXkhMRGR3DrwBS+Wbmdxdl7/I4jIj46UtF9ObAF+MbMXjazMwEtBScinlqycS9fL9/GzwckExetXm6p+W7sl0iDOpE8MznD7ygi4qPDFt3OuY+dcz8DOgPpwK+B5mb2opmdU035RCTEPP31KurHRDC8f5LfUUSqRL2YSG45PZmvl29lyca9fscREZ9U5kbKXOfcO865C4HWwELgPs+TiUjIWbhhD18v38bIQSk0qBPpdxyRKjO8fxL1YyJ4+muN7RYJVce0prJzbpdzbpRzrlJLwJvZUDNbaWYZZvaTQt3MBpnZAjMrMrMryh0rNrOFgcf4Y8kpIjXTk5NW0Tg2ipsGJPsdRaRKNagTychBKXy9fCsLN2hst0goOqai+1iYWTjwPHAekApcY2ap5ZqtB24C3q3gJfKcc90Dj4u9yikiwWHu2l1MW7Wd2wenaCy31Eo3DUimcWwU//5qpd9RRMQHnhXdQB8gwzmX6ZwrAMYCl5Rt4Jxb65xbDJR4mENEaoB/f7WS+HrR3HBakt9RRDwRFx3BHYPbMX31DuZk7vQ7johUMy+L7gRK5/Y+JDuwr7JizGyemc02s0urNpqIBJOZGTuYnbmLO9PaUSdKq09K7XX9aYk0qxfNv79ahXPO7zgiUo28/B1uRdMLHss3TFvn3CYzSwGmmNkPzrk1P3oDs5HASID4+HjS09OPO6zUTjk5ObougpxzjkfmHKRxjJFwcC3p6es8f09dF1KR6rouzmnteHv5Lp7/aDJdm2ooVbDT94VUFS//b88G2pTZbg1squzJzrlNgT8zzSwd6AGsKddmNDAaoFOnTi4tLe3EEkutk56ejq6L4PbNym1k7JnLI5d15Zy+1bMYjq4LqUh1XRf9ior55l9TmbQlirsuH4CZlsAIZvq+kKri5fCSuUAHM0s2syjgaqBSs5CYWSMziw48bwoMAJZ5llREfOGc46lJq2jdqA5X9mpz9BNEaoHoiHB+eWYHFmXvZdKyrX7HEZFq4lnR7ZwrAu4GJgLLgXHOuaVm9pCZXQxgZqeaWTZwJTDKzJYGTu8CzDOzRcA3wGPOORXdIrXMF0u2sDh7L788swNREV72AYgEl2E9E0hqUpcnJ62ipERju0VCgaeDyZxzE4AJ5fY9UOb5XEqHnZQ/bybQzctsIuKvwuISnpi4kk7N6zGs50++BkRqtYjwMH59dkd+OXYh4xdt4tIexzLPgIjUROpaEhFfvD93A1k7cvnDeZ0ID9OYVgk9F53citSW9fnXVyvJLyr2O46IeExFt4hUu9z8Ip7+ejV9khszpFMzv+OI+CIszPjj+Z3J3p3HW7O8n7VHRPyloltEqt2YGVnsyMnnvvM6a+YGCWkDO8QzsENTnvsmg715hX7HEREPqegWkWq1MyefUdMyGXpSC3q2beR3HBHf/WFoZ/YcKOTF9DVHbywiNZaKbhGpVs99k0FeYTG/G9rJ7ygiQaFrQgMu7d6K177NYtOePL/jiIhHVHSLSLVZv/MAb89ex1W929AuPs7vOCJB4zfndMI5eGrSKr+jiIhHVHSLSLX511crCQ8zfnVWB7+jiASVNo3rcmO/RD5akM3KLfv9jiMiHlDRLSLVYv663YxftIkRA1NoXj/G7zgiQeeuIe2JjY7gsS+W+x1FRDygoltEPFdS4njo82U0rx/N7YPb+R1HJCg1io3i7iHt+Wbldqau2u53HBGpYiq6RcRznyzcyKINe/jD0M7ERnu6EK5IjXbTgCQSm9Tl4c+XUVhc4nccEalCKrpFxFMHCor455crOKV1Ay7trqWuRY4kOiKcv1yQSsa2HN6erQVzRGoTFd0i4qmX0tewdV8+D1yUSpiWexc5qrO6NGNgh6Y8NWkVu3IL/I4jIlVERbeIeGbjnjxGTcvkolNa0Suxsd9xRGoEM+P+C1PJLSjWFIIitYiKbhHxzD+/WAHAfed19jmJSM3SsXk9ru/blnfmrGPFln1+xxGRKqCiW0Q8MX/dLsYv2sRtg1JIaFjH7zgiNc6vz+5I/TqRPPTZMpxzfscRkROkoltEqlxRcQl/+WQpLerHcJumCBQ5Lg3rRvHrszoyc81OJi7d6nccETlBKrpFpMq9OWsdyzfv48GLUjVFoMgJuK5vWzo1r8fDny/jQEGR33FE5ASo6BaRKrVt30GenLSKQR3jGdq1hd9xRGq0iPAwHr60Kxv35PHM5Ay/44jICVDRLSJV6u//W05BUQl/u/gkzDRFoMiJ6pPcmCt6teaV6Zms2rrf7zgicpxUdItIlZmZsYPxizZxe1o7kpvG+h1HpNb443mlq7n+5ZMluqlSpIZS0S0iVaKgqIT7P11Cm8Z1uDNNN0+KVKUmcdHcd15nvsvaxUcLNvodR0SOg4puEakSo6etYc32XB66uCsxkeF+xxGpdX7Wuw092zbkHxOWs+eAVqoUqWk8LbrNbKiZrTSzDDO7r4Ljg8xsgZkVmdkV5Y4NN7PVgcdwL3OKyInJ2JbDM5MzOL9bC4Z0buZ3HJFaKSzMeOSybuzNK+TRCSv8jiMix8izotvMwoHngfOAVOAaM0st12w9cBPwbrlzGwMPAn2BPsCDZtbIq6wicvxKShz3fbSYOlHh/O3irn7HEanVurSsz60Dk3l/3ga+zdjhdxwROQZe9nT3ATKcc5nOuQJgLHBJ2QbOubXOucVASblzzwUmOed2Oed2A5OAoR5mFZHj9Pacdcxbt5u/XNCF+HrRfscRqfV+fVZHkpvGct9/F2vubpEaxMuiOwHYUGY7O7DP63NFpJps3JPHP79YwcAOTbmiV2u/44iEhJjIcB4b1o0Nu/J4YuJKv+OISCV5uVRcRRP0Vnaeo0qda2YjgZEA8fHxpKenVzqchIacnBxdFx5xzvHUgnwKi4u5uGUuU6dO9TtSpem6kIrUtOvijLYRvP7tWloVbaFDI9287JWadl1I8PKy6M4G2pTZbg1sOoZz08qdm16+kXNuNDAaoFOnTi4tLa18Ewlx6enp6LrwxsffZ7N4+yLuvzCVK09P9jvOMdF1IRWpaddF735FnPvUNMZmhvG/ewZq1iCP1LTrQoKXl8NL5gIdzCzZzKKAq4HxlTx3InCOmTUK3EB5TmCfiASBzXvzeODTpfRs25Cb+if5HUckJMVFR/CPYd1Ysz2X/0xe7XccETkKz4pu51wRcDelxfJyYJxzbqmZPWRmFwOY2almlg1cCYwys6WBc3cBD1NauM8FHgrsExGflZQ4fvfBYopLHE9e1Z3wMC31LuKXwR3juap3a0ZNXcO8tfoxKRLMvBxegnNuAjCh3L4HyjyfS+nQkYrOHQOM8TKfiBy7t+esY0bGDh65rCtJWupdxHcPXHQSszJ3cu+4RUz45UDioj390S4ix0krUopIpa3ZnsM/JixncMd4ru3T1u84IkLpMJMnr+pO9u4DPPzZMr/jiMhhqOgWkUopKi7h3nGLiI4I5/ErTsZMw0pEgsWpSY25fXA73p+3gYlLt/gdR0QqoKJbRCrl2SkZLNqwh79f2pXm9WP8jiMi5fzqrI6ktqzPH//7A9v35/sdR0TKUdEtIkc1a81Onp2ymmE9ErjolFZ+xxGRCkRFhPH01d3JyS/itx8soqSksktjiEh1UNEtIke0MyefX479nqQmsTx8aVe/44jIEXRsXo/7L+jC1FXbGTUt0+84IlKGim4ROaySEsdvPljEnrxCnru2J7GaFUEk6F1/WiIXdGvJv75aqWkERYKIim4ROaxXZmSSvnI791/QhdRW9f2OIyKVYGY8enk3Wjeqwy/e+57duQV+RxIRVHSLyGF8v343j3+5kvO6tuD60xL9jiMix6B+TCTPX9uTnTkF/Ebju0WCgopuEfmJ7fvzuePtBbRsGMNjl2t6QJGaqGtCA/5yYRemrNim8d0iQUBFt4j8SGFxCXe9s4A9eQWMur43DepE+h1JRI7TDaclcsHJLXli4gqmrtrudxyRkKaiW0R+5JH/Lee7tbv45+Unaxy3SA1nZjxxxcl0bF6PX7y7gLU7cv2OJBKyVHSLyP/5cH42r89cy4iByVzSPcHvOCJSBepGRfDyjb0JCzNGvDmPnPwivyOJhCQV3SICwOLsPfzp4x/o364Jfxja2e84IlKF2jSuy/PX9iRzRy73vr9QN1aK+EBFt4iwcU8et7wxj/i4aJ67ticR4fpqEKltBrRvyp/O78JXy7by9OTVfscRCTla6UIkxO07WMjNr33HwcJi3r21L41jo/yOJCIe+fmAJFZs3sczk1fTtnFdrujV2u9IIiFDRbdICCssLuHOtxeQuT2XN37ehw7N6/kdSUQ8ZGY8clk3Nu89yH0fLaZF/RhO79DU71giIUG/QxYJUc45/vLxEmZk7ODRYd0Y0F4/eEVCQVREGC9c35P2zeK4/e35LN+8z+9IIiFBRbdIiHp2Sgbvz9vAL85oz5W92/gdR0SqUf2YSF67+VTioiO4+bW5bN6b53ckkVpPRbdICHr92yyenLSKYT0SuPfsjn7HEREftGxQhzE3nUpOfhHDx3zHrtwCvyOJ1GoqukVCzIfzs/nrZ8s4J7U5j1+hJd5FQllqq/qMvrEX63Ye4MYxc9h3sNDvSCK1lopukRDy5ZLN/P7DRZzevinPXttDUwOKCP3bNeWl63uxcst+bn5tLgcKtHiOiBf0E1ckRExbtZ173ltI9zYNGXVDL6Ijwv2OJCJBYkjnZvzn6h58v343I96cx8HCYr8jidQ6nhbdZjbUzFaaWYaZ3VfB8Wgzez9wfI6ZJQX2J5lZnpktDDxe8jKnSG03ZcVWbn1zHu2axfHaTX2IjdZsoSLyY+d3a8kTV5zCtxk7uf3t+Sq8RaqYZ0W3mYUDzwPnAanANWaWWq7ZLcBu51x74Cngn2WOrXHOdQ88bvcqp0ht9+WSLdz21nw6Na/HeyP60qBupN+RRCRIXd6rNY8O68bUVdu55Q0NNRGpSl72dPcBMpxzmc65AmAscEm5NpcAbwSefwicabqrS6TKfL54E3e9u4CuCQ14+9a+NKyr1SZF5Miu6dOWf195CrPW7OTGV7/TzZUiVcScc968sNkVwFDn3K2B7RuAvs65u8u0WRJokx3YXgP0BeKApcAqYB/wF+fc9AreYyQwEiA+Pr7XuHHjPPksUnPl5OQQFxfndwxfzNhYyKs/FNChURi/7hVDnQj9e/aQUL4u5PB0XfzYd1uKGLUon7b1wvhN7xjiokLzO0TXhVRkyJAh851zvY/lHC8Hdlb0f2f5Cv9wbTYDbZ1zO82sF/CJmZ3knPvRslnOudHAaIBOnTq5tLS0E08ttUp6ejqhdl0453ghfQ2v/LCS/u2a8PKNvTWGu5xQvC7k6HRd/Fga0OuUrdzxzgL+sySM124+ldaN6vodq9rpupCq4uXwkmyg7DJ3rYFNh2tjZhFAA2CXcy7fObcTwDk3H1gDaAUPkaMoLnHc/+kSnpi4kku6t+L1m3XTpIgcvzO7NOf1m09ly76DXPbCTJZs3Ot3JJEay8uiey7QwcySzSwKuBoYX67NeGB44PkVwBTnnDOz+MCNmJhZCtAByPQwq0iNl1dQzO1vz+ft2eu5fXA7nrqqO1ERmhVURE5M/3ZN+eiO/kSFh3HVqFl8s2Kb35FEaiTPfiI754qAu4GJwHJgnHNuqZk9ZGYXB5q9CjQxswzgXuDQtIKDgMVmtojSGyxvd87t8iqrSE23aU8eV42axdfLt/LQJSdx33mdCQsLzfGXIlL1Ojavx8d39iclPpZb35zHW7PX4dU9YSK1lae/d3bOTQAmlNv3QJnnB4ErKzjvI+AjL7OJ1BZzMndy5zsLyC8q4eUbenNWanO/I4lILdSsfgzvj+zHL977nvs/WcKS7L387ZKTiInUQlsilaHfPYvUUM453py1lutemUODOpF8ctcAFdwi4qnY6AhevrE3dw9pz/vzNnDVqFls3JPndyyRGkFFt0gNlJNfxG8+WMQDny5lUMd4Prl7AO2baUorEfFeeJjx23M7MeqGXmRuz+WiZ2cwM2OH37FEgp6KbpEa5ofsvVz4zHQ++X4jvzyzA6/c2Jv6MVplUkSq17knteDTuwfQODaK616dwz+/XEFBUYnfsUSClopukRqipMTx8rRMhr34LflFJbw34jR+fXZH3TApIr5pFx/H+LsH8LPebXgxfQ2XvziTzO05fscSCUoqukVqgPU7D3DDmDk8MmE5Z3Ruxhe/HEjflCZ+xxIRoW5UBI9dfjIvXd+TDbsPcMEzM3hnjmY3ESlPq2aIBLHiEsfrM9fyr4krCQ8zHh3WjatPbYOZerdFJLgM7dqS7m0a8ZsPFvLnj5cwfuEmHh3WjZR43W8iAurpFglaq7bu58qXZvLw58s4LaUxX/16ENf0aauCW0SCVosGMbz18748NqwbyzbvY+h/pvP8NxkUFmust4h6ukWCzN68Qp7+ehVvzlpHvZgInvrZKVzaPUHFtojUCGFhxtV92nJG52Y8OH4pT0xcyfiFm3jgolQGtG/qdzwR36joFgkSxSWOcfM28MTElew+UMA1fdrym7M70iQu2u9oIiLHrFn9GF68vhcTl27h4c+Xcd0rczg7tTl/Or8LyU1j/Y4nUu1UdIv4zDnH18u38e+vVrJiy35OTWrEgxf1oWtCA7+jiYicsHNPasHgjvG8OiOLF77J4JynpnJjvyTuSGtHU3UqSAhR0S3io28zdvDExJUs3LCHpCZ1eeaaHlx0cksNJRGRWiUmMpy7hrTnyl6t+ddXK3nt2yze+249w/snMXJgCo1io/yOKOI5Fd0i1aykxJG+ahsvpWfy3dpdtGoQw2PDunF5r9ZEhuveZhGpvZrVj+HxK05h5KB2/Gfyal6auoa3Zq3jpv5JDO+fRHw99XxL7aWiW6SaFBSVMH7RJkZPW8OqrTkkNKzDXy9K5Zq+bYmOCPc7nohItWnfLI5nr+nB3UPa85/Jq3g+PYPR0zO5vGcCt5yeQvtmmmZQah8V3SIe27gnj/fmrOf9eRvYvj+fzi3q8dTPTuHCk1upZ1tEQlqnFvV44bpeZG7P4dUZWXw4P5v/197dx8hR33ccf3939/bu1vds+872ncEGP4DB2IRHF5LaPIUoSZ2oRZC2FFVIqBJPSaNWoVKbqBISUaO2VEGRKNAGGoEoQYmDUIBgGxoVsCE4Bhs/YYN99tln+3z2ne9xd7/9Y8b23nl9ZxfGs7f+vKTRzO83v5n9ru/n3e/O/Gbm2TW7WDZ/Kt8K74CS0ueklAkl3SIRGMrmeWPLfp5bs5NVmztx4Ib5zdy55Hz+cN5UjdkWESlwwdQaHv7mQv765nk8/danPLtmJ/c88x7NtZXcdmUbd1x1HjObMnGHKfKZKOkW+Zzk8867nx7iF+t28/IHHXT3DTO1tpJ7l83h9qtm0taoLwwRkbFMrqnkOzfP4/4b5rByUyfPrd3FT1Z/zGOrPubq2U18/bLpfGXhdN31RCYkJd0in8FwLs/aT7r4zcZOXtmwl93d/VRXJLnlkha+sbiV6+dO0RASEZEzlEomuOWSadxyyTT2dPfz3++286v1e/j7X27g+ys2sOTCyXztshnceHEzzbVVcYcrclqUdIucof09g/zvxwdYuamTVZs6OTKQJZ1KcP2cKfzNl+dz84IWJlXqv5aIyOdhRkM1D940lwdunMPmfT289PsOXlq/h4de/ACAha31LJs/lWUXNXNZWwPJhIbvSWlSZiAyju6+Id7e3sVbHx/gre0H2bKvF4DJk9J8+ZJp3LSghS/OnUImrf9OIiJRMTMumlbHRdPq+O4t8/ioo4dVm4ODHz9etY1/W7mNhkwFV81q4prZTVwzezILZtQpCZeSoSxBpMBwLs/mvT2s29XN73d1s25XN9v29+IO1RVJrpzVyDcvb2PJhZNZ2FqvD3MRkRiYGQtm1LFgRh33LptDd98Qb2zZz2+3HmDNJ128tnEfALWVKa6Y1ciitgYua6tnYWs9zXUajiLxUNIt5yR3Z3/vIFv29rJ5Xw9b9vawaV8PmzqOMJjNA8GR7EUzG/j6ohksuXAyi9oaSKc0PltEpNQ0ZNIsX9zK8sWtAOw9PMA7Ow6yZkcXa3Z08caW/bgHbVvqKlnYWpx16D4AAAtOSURBVM+C6XXMaallztQaLpg6iaoKPS9BoqWkW8qWu9M9mOe9T7vY2dXHzoP9fNp1lF1dfWzr7OVQ3/DxtlNq0sxrqeXOa89n8XkNLGproK2xWrf2ExGZgKbVV41Iwo8OZtnYcYT17Yf5cPdh1rd3s3JTJ/kwEU8YzGzKMLe5htlTJtHWmGFmUzUzGzMMZj3GdyLlJNKk28xuBR4FksAT7v7IqPWVwNPAFcBB4HZ3/yRc9xBwN5ADHnD3V6KMVSaOfN45MjDMob5h9vcM0tkzwL4jg3QeGWDfkWB5X88AHd0D9A/nYNVbAJjB9Loq2poy3HrpdOa31DBvWi3zWmp1+ykRkTI2qTLFVbOauGpW0/G6geEcnxw8ytZ9vWzrDKatnT38z9YDx894HjPl7ddoa8zQUldJc20VU2sraa6tDOdBeXJNWnerkjFFlnSbWRJ4DLgZaAfWmtkKd99Y0Oxu4JC7zzGzO4AfAreb2QLgDuASYAbwGzOb5+65qOKV6Lk7g9k8fUM5jg5mg/lQlr7BHL2DWfqGshwdytE3GMx7Bobp7humu2+IQ+G8u3+Yw/3Dx08TFkqnErTUVdJSW8XF0+pYNr+ZwYO7ufHaRZzXlKG1oVqnD0VEBICqiuTxCzMLHRt+uKurn/ZDfbz53gYqGlpoP9TP9v1HeWdHF90FZ0oLTUonqa+uoD6Tpr46RX11BQ3VaeozFdRXV1BTmSKTTpJJp8hUJpmUPlY+UZepSOopnGUqyiPdVwPb3H07gJk9BywHCpPu5cAPwuUXgB9bcD5/OfCcuw8CO8xsW7i/t071Yn1Z59UNeynMxU4kZn5SXbF2Pm67kzO9MbctaO7H1/lJdYxo52Nse3rtKPIa7pB3J5f3cH6ifKLOybmTz49a704uV7AurM/mnKFcnqFsMA0eX84xlMsznB25fig38sjBeGoqUzRkKmjIVNCYSdPWWE1jJk1jJvhAa6iuYGptJS11VbTUVVJfXXHScJDVqztZOr/5jF5XRETOXWZGc20VzbVVXHF+I/XdW1m69LIRbQazOQ70DtF5ZCA84zpI19EhDvcHB4sO9w9zuH+IHQeO0t3XTXf/MEPZ0/8OTCWMimSCdCqckgkqw+Xj9YXrUwmSZiQT4WRGImEkE4T1CZIJgrqi7YxUwo5/hxrB2WEL/z2OLWNWsM4K2gRlCrc5Vl/QltHrCupP+fcY9+/1/9/6s772mYoy6W4FdhWU24FrTtXG3bNmdhiYHNa/PWrb1tEvYGb3APcApKfN4Z5n3vvcgj/XGMGYNrNgngjLJyYbWQaSifCDIQGpBKTMqEzApARUpCBVXbAukaQikaQiAZUpoypZME8aVSmoShqV4TydDF4zkAP6wymUBXog3wMde6DjFO+rt7eX1atXR/bvJhOT+oUUo34hxYzXL9JAG9CWBGrCaYQUkGIo5wxkYTDnDOSC+WBBuXD9cB6G804272TzObL5XFiGbBb6B5yePAznIRu2DQ6uBVMunLs7eUbXBcty9kWZdBf7gTD6z3yqNqezLe7+OPA4wOw58/1X918f7LRgawt3NaLORq4bWXdyO8ZtZ0XqTv0aFGtXsNKKtC/6Pk4Ob8x2iYJftMGvW47/4i3XCwZXr17N0qVL4w5DSoz6hRSjfiHFlGu/KDxznfcgwXcHPDib7h4kXu4ezsOz7Mfri7QrONN/qn0won7sGP3k1G/k+jFWR7lvgIU/HHt9MVEm3e3AzIJyG7DnFG3azSwF1ANdp7ntCOkkXNpa/1ljFhERESl7iYSRwNClTmdPlCP11wJzzWy2maUJLoxcMarNCuCucPlPgJUeDHpeAdxhZpVmNhuYC6yJMFYRERERkchEdqQ7HKN9H/AKwS0Dn3L3DWb2j8C77r4CeBJ4JrxQsosgMSds9zzBRZdZ4F7duUREREREJqpI79Pt7i8DL4+q+4eC5QHgtlNs+zDwcJTxiYiIiIicDboRpIiIiIhIxJR0i4iIiIhETEm3iIiIiEjElHSLiIiIiERMSbeIiIiISMSUdIuIiIiIRExJt4iIiIhIxMzHe7j8BGFmPcDmuOOQkjMFOBB3EFJy1C+kGPULKUb9QoqZ7+61Z7JBpA/HOcs2u/uVcQchpcXM3lW/kNHUL6QY9QspRv1CijGzd890Gw0vERERERGJmJJuEREREZGIlVPS/XjcAUhJUr+QYtQvpBj1CylG/UKKOeN+UTYXUoqIiIiIlKpyOtItIiIiIlKSyiLpNrNbzWyzmW0zs+/FHY/Ez8xmmtkqM/vIzDaY2YNxxySlw8ySZva+mb0UdyxSGsyswcxeMLNN4efGkrhjkviZ2XfC75APzexZM6uKOyY5+8zsKTPrNLMPC+qazOw1M9sazhvH28+ET7rNLAk8BnwFWAB8y8wWxBuVlIAs8F13vxi4FrhX/UIKPAh8FHcQUlIeBX7t7hcBi1D/OOeZWSvwAHClu18KJIE74o1KYvKfwK2j6r4HvO7uc4HXw/KYJnzSDVwNbHP37e4+BDwHLI85JomZu3e4++/C5R6CL9DWeKOSUmBmbcBXgSfijkVKg5nVAV8CngRw9yF37443KikRKaDazFJABtgTczwSA3d/E+gaVb0c+Gm4/FPgG+PtpxyS7lZgV0G5HSVXUsDMZgGXA+/EG4mUiH8F/hbIxx2IlIwLgP3Af4TDjp4ws0lxByXxcvfdwI+AnUAHcNjdX403KikhLe7eAcGBPqB5vA3KIem2InW6JYsAYGY1wM+Bb7v7kbjjkXiZ2deATnd/L+5YpKSkgC8AP3H3y4GjnMapYilv4Rjd5cBsYAYwycz+PN6oZCIrh6S7HZhZUG5Dp38EMLMKgoT7Z+7+YtzxSEm4DvgjM/uEYCjaDWb2X/GGJCWgHWh392Nnw14gSMLl3HYTsMPd97v7MPAi8AcxxySlY5+ZTQcI553jbVAOSfdaYK6ZzTazNMFFDitijkliZmZGMD7zI3f/57jjkdLg7g+5e5u7zyL4rFjp7jpydY5z973ALjObH1bdCGyMMSQpDTuBa80sE36n3IgusJUTVgB3hct3Ab8cb4NUpOGcBe6eNbP7gFcIrix+yt03xByWxO864E7gAzNbF9b9nbu/HGNMIlK67gd+Fh682Q78ZczxSMzc/R0zewH4HcEdsd5HT6c8J5nZs8BSYIqZtQPfBx4Bnjezuwl+oN027n70REoRERERkWiVw/ASEREREZGSpqRbRERERCRiSrpFRERERCKmpFtEREREJGJKukVEREREIjbhbxkoInIuMrPJwOthcRqQI3iUOUCfu+shHiIiJUS3DBQRmeDM7AdAr7v/KO5YRESkOA0vEREpM2bWG86XmtkbZva8mW0xs0fM7M/MbI2ZfWBmF4btpprZz81sbThdN87+p5vZm2a2zsw+NLMvno33JSIykWl4iYhIeVsEXAx0ETxp8Ql3v9rMHiR4CuO3gUeBf3H335rZeQRP+L14jH3+KfCKuz9sZkkgE+k7EBEpA0q6RUTK21p37wAws4+BV8P6D4Bl4fJNwAIzO7ZNnZnVunvPqfYJPGVmFcAv3H1dNKGLiJQPDS8RESlvgwXL+YJynhMHXhLAEndfHE6tYyTcuPubwJeA3cAzZvYXEcQtIlJWlHSLiMirwH3HCma2OJxfbWZPj25sZucDne7+78CTwBfOVqAiIhOVkm4REXkAuNLM1pvZRuCvwvrzgP4i7ZcC68zsfeCPCcaEi4jIGHTLQBERKcrM/gl4xt3Xxx2LiMhEp6RbRERERCRiGl4iIiIiIhIxJd0iIiIiIhFT0i0iIiIiEjEl3SIiIiIiEVPSLSIiIiISMSXdIiIiIiIRU9ItIiIiIhKx/wNYX7ej7H0k7AAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 864x360 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Define figure size\n",
    "rcParams['figure.figsize'] = 12, 5\n",
    "\n",
    "# Initial parameters\n",
    "tmax=10.0                     # maximum time\n",
    "nt=1000                       # number of time sample\n",
    "a=1                           # half-width      \n",
    "dt=tmax/nt                    # defining dt\n",
    "t0 = tmax/2                   # defining time shift t0\n",
    "\n",
    "time=np.linspace(0,tmax,nt)   # defining time\n",
    "\n",
    "# Define gaussian function          \n",
    "f=(1./np.sqrt(2*np.pi*a))*np.exp(-(((time-t0)**2)/(2*a)))\n",
    "    \n",
    "# Plotting of gaussian\n",
    "plt.plot(time, f)\n",
    "plt.title('Gaussian function')\n",
    "plt.xlabel('Time, s')\n",
    "plt.ylabel('Amplitude')\n",
    "plt.xlim((0, tmax))\n",
    "plt.grid()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "Next, we calculate the numerical derivative using finite-differences with the forward operator:\n",
    "\n",
    "\\begin{equation}\n",
    "\\biggl(\\frac{\\partial f(t)}{\\partial t}\\biggr)^+ \\approx \\frac{f(t+dt)-f(t)}{dt}. \\nonumber\n",
    "\\end{equation}\n",
    "\n",
    "backward operator:\n",
    "\n",
    "\\begin{equation}\n",
    "\\biggl(\\frac{\\partial f(t)}{\\partial t}\\biggr)^- \\approx \\frac{f(t)-f(t-dt)}{dt}, \\nonumber\n",
    "\\end{equation}\n",
    "\n",
    "and central operator:\n",
    "\n",
    "\\begin{equation}\n",
    "\\biggl(\\frac{\\partial f(t)}{\\partial t}\\biggr) \\approx \\frac{f(t+dt)-f(t-dt)}{2dt}. \\nonumber\n",
    "\\end{equation}\n",
    "\n",
    "To test the accuarcy of these approaches, we compare them with the analytical derivative:\n",
    "\n",
    "\\begin{equation} \n",
    "\\frac{\\partial f(t)}{\\partial t}=-\\dfrac{t-t_0}{a}\\dfrac{1}{\\sqrt{2\\pi a}}e^{-\\dfrac{(t-t_0)^2}{2a}} \\nonumber\n",
    "\\end{equation} "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "code_folding": []
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAtwAAAFNCAYAAAAze7gSAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJzs3XdcleX/x/HXxREFFfceiVqKgym4cpKaOcmRM0fD3c4yy5HasmEafiut1JQQd2VaTnKnIGgqSpY4caEiyOZcvz/A80MFwXG4GZ/n43EennOP63rf3B74nPtc930rrTVCCCGEEEII67AxOoAQQgghhBAFmRTcQgghhBBCWJEU3EIIIYQQQliRFNxCCCGEEEJYkRTcQgghhBBCWJEU3EIIIYQQQliRFNxCCJGLlFKPKKVilVKmXOpvqlJqyQOsv14pNfRhZsrQdqxSqo412hZCiLxECm4hhLACpVSEUio+vai8+aimtT6ltS6ptU69jzaHKaV2WCNvVrTWT2mtFz1oO0qpQKXUC7e1XVJr/d+Dti2EEHldEaMDCCFEAdZda70ppwsrpRSgtNZmK2bKd1mEECK/kyPcQgiRi5RSjkoprZQqkv46UCn1gVJqJxAH1Ek/kv2fUipGKXVCKTVIKdUA+AZokX60/FoW7ddWSv2Zvu5GoMJt85srpXYppa4ppQ4opdplmJdZlkCl1AtKqWLp6zTOsHzF9KP4lZRSZZVSa5VSl5RSV9Of10hf7gOgNeCbnt03fbpWSj2anul8xmE2SqmnlVIH05/bKKUmKKX+VUpFKaWWKaXKPfjeEEKI3CEFtxBCGO9ZYATgAFwC5gBPaa0dgJZAqNY6DBgF7E4filEmi7Z+AoJJK7SnA5bx10qp6sBvwAygHPAmsFIpVTGLLCdvTtRaJwKrgAEZln0G+FNrfZG0vycLgFrAI0A84Ju+7rvAdmBcevZxGQNrrfcANwDvDJMHpm8LwMuAD9AWqAZcBeZmsf1CCJHnSMEthBDWsyb9qPA1pdSauyy3UGt9WGudAqQAZqCxUspeax2ptT6ck86UUo8AXsAkrXWi1nob8GuGRQYD67TW67TWZq31RiAI6JJZFq118m1d/MStBbelKNZaR2mtV2qt47TWMcAHpBXIOeV/s22llEN6Jv/0eSOBd7XWZ9IL/6lAn5vfEgghRF4nBbcQQliPj9a6TPrD5y7Lnb75RGt9A+hH2tHsSKXUb0oppxz2Vw24mt7GTSczPK8F9M3wIeAa0AqomlmWTGwB7JVSzZRStQA3YDWAUqq4UupbpdRJpdR1YBtQ5h6uxvIT0EspVQzoBezXWt/MXgtYnSFzGJAKVM5h20IIYSgpuIUQwnj6lhda/6G17khaIXwUmJ/ZcpmIBMoqpUpkmPZIhuengcUZPgSU0VqX0Fp/nFWW23KZgWWkHYkeCKxNP5oN8AZQH2imtS4FtEmfrnKSXWt9hLQPB09x63CSm7mfui23ndb67N3aFEKIvEIKbiGEyEOUUpWVUj3Si+ZEIJa0o7kAF4AaSqmima2bfkQ4CHhfKVVUKdUK6J5hkSVAd6XUk0opk1LKTinV7ubJjTn0E2lH4Adxa1HsQNq47WvpJzROuW29C0B219z+ibTx2m2A5RmmfwN8kH5U/ebJmj3vIbMQQhhKCm4hhMhbbEg7WnwOuELaOOgx6fO2AIeB80qpy1msPxBolr7uFODHmzO01qeBnsBE0k7OPA2M5x7+Fmit/yLtBMdqwPoMs74E7IHLwB7g99tWnU3auOurSqk5WTTvD7QDtmitM27fbOAXYINSKia9/WY5zSyEEEZTWmf3DaUQQgghhBDifskRbiGEEEIIIaxICm4hhBBCCCGsSApuIYQQQgghrEgKbiGEEEIIIaxICm4hhBBCCCGsqMDdFrdChQra0dHRkL5v3LhBiRIlsl9Q5Guynws+2ceFg+znwkH2c+Fg1H4ODg6+rLWumN1yBa7gdnR0JCgoyJC+AwMDadeunSF9i9wj+7ngk31cOMh+LhxkPxcORu1npdTJnCwnQ0qEEEIIIYSwIim4hRBCCCGEsCIpuIUQQgghhLCiAjeGWwghhBDiJqUUJ06cICEhwegowopKly5NWFiY1dq3s7OjRo0a2Nra3tf6UnALIYQQosAqUaIEDg4OODo6opQyOo6wkpiYGBwcHKzSttaaqKgozpw5Q+3ate+rDRlSIoQQQogCy2QyUb58eSm2xX1TSlG+fPkH+pZECm4hhBBCFGhSbIsH9aD/h6TgFkIIIYSwIpPJhJubm+URERFBYGAgpUuXxt3dnfr169OmTRvWrl2b6fqJiYl06NABNzc3AgICcjn9ndq1a2fYPU8yioiI4KeffjI6Ro7IGG4hhBBCCCuyt7cnNDT0lmkRERG0bt3aUmSHhobi4+ODvb09TzzxxC3LhoSEkJycfEcbd5OamorJZHrg7CkpKRQpYly5eLf+bxbcAwcOfCjtWZMc4RZCFGqpqQn8vvlbPh7XjxeHePPoqz0YtGIGj77Uh/e7tmP6G+/w555wzGajkwohCjI3NzcmT56Mr6/vLdMvXrzI4MGDCQ0Nxc3NjX///ZfNmzfj7u6Os7Mzzz33HImJiUDa3banTZtGq1at8PPzo0mTJgAcOHAApRSnTp0CoG7dusTFxfHrr7/SrFkz3N3d6dChAxcuXABg6tSpjBgxgk6dOjFkyBDi4+Pp378/Li4u9OvXj/j4+Ey34W653n77bZo2bUrTpk05fvw4AJcuXaJ37954eXnh5eXFzp07M+3/5ocTDw8PPDw82LVrFwATJkxg+/btuLm54evrS0JCAsOHD8fZ2Rl3d3e2bt0KwMKFC+nbty/du3enU6dOD22f3ROtdYF6NGnSRBtl69athvUtco/s5/wvJSVBr145Tc9qV1fvrlZEpyi0Br2hDpqpaY9SE9KmadAxtuiAumX1hKf76ZCDV4yOLx4SeS8XDvv37zc6graxsdGurq7a1dVV+/j4aK3T/v917dr1luVCQkK0k5PTHetnXDY+Pl7XqFFDHzt2TGut9bPPPqtnzZqltda6Vq1a+pNPPrGs17BhQx0dHa2/+uor7enpqZcsWaIjIiJ08+bNtdZaX7lyRZvNZq211vPnz9evv/661lrrKVOmaA8PDx0XF6e11vrzzz/Xw4cP11prfeDAAW0ymfS+fftuyZhdrhkzZmittV60aJFlWwYMGKC3b9+utdb65MmTlm2/vf8bN27o+Ph4rbXW4eHh+matl/Hncv36df3ZZ5/pYcOGaa21DgsL0zVr1tTx8fF6wYIFunr16joqKirrnZQDR44cuWMaEKRzUJ/KkBIhRKGRkhLP4q/HUvJ/K3j6WAxFdNr0JBsIrViKyLK16G7zJCnxcZRW8Wxw2kCd85d59Foiz/x7lWf+DeDYn8t5x7UTz/5vOQ2dShq7QUKIexIYaJ2TJ9u103edn9mQksyk1W93d+zYMWrXrk29evUAGDp0KHPnzuXVV18FoF+/fpZlW7Zsyc6dO9m2bRsTJ07k999/R2tN69atAThz5gz9+vUjMjKSpKSkWy5516NHD+zt7QHYtm0bL7/8MgAuLi64uLjcc64BAwZY/n3ttdcA2LRpE0eOHLG0cf36dWJiYu7oPzk5mXHjxhEaGorJZCI8PDzTn82OHTt46aWXAHBycqJWrVqWZTt27Ei5cuWy/flaixTcQogCT2tN4N/zGeQ3m95/H+Gro5BsA7/XeoQEnxd46r1XaVLegSbAECAwMJB27drBh2nrRx3cwq4Zk3D5Yy/1r6QwPfB3Gn/kRNM6s/nh3d4YOLxRCFGAhISE0KBBg7suk11RXqJECcvz1q1bs337dk6ePEnPnj355JNPUErRrVs3AF566SVef/11evToQWBgIFOnTs20Hcj+Kh3Z5cq4/s3nZrOZ3bt3WwrrrLZj1qxZVK5cmQMHDmA2m7Gzs7vnDLdvT26TPxNCiAItLv4sYz7owyKb/VA8ie8alMYzsj5NPv2ezh0a56iN8i7edF/mjTkpjs0TnmHnvu0cq3OWY/Rh0wtj2TntC2o/UtTKWyKEeFDZHYk20sGDB5k+fTrffffdXZdzcnIiIiKC48eP8+ijj7J48WLatm2b6bJt2rThvffeo02bNtjY2FCuXDnWrVvHRx99BEB0dDTVq1cHYNGiRVn22aZNG/z8/Gjfvj2HDh3i4MGD95wrICCACRMmEBAQQIsWLQDo1KkTvr6+jB8/HsAyTv120dHR1KhRAxsbGxYtWkRqaioADg4OliPiGXN6e3sTHh7OqVOnqF+/Pvv377/rzzQ3yEmTQogCK2T/PP7wcGLmrD1UiUui7Ml+rOxxhKEhf9E4h8V2RjZFi/PEF2t5Y8NhepZsD2YTrilz2e1dg42bI62wBUKIgmz79u2WywKOHTuWOXPm3HGFktvZ2dmxYMEC+vbti7OzMzY2NowaNSrTZR0dHYG0QhSgVatWlClThrJlywJpJyf27duX1q1bU6FChSz7HD16NLGxsbi4uDBz5kyaNm16z7kSExNp1qwZs2fPZtasWQDMmTOHoKAgXFxcaNiwId98802m/Y8ZM4ZFixbRvHlzwsPDLUerXVxcKFKkCK6urvj6+jJmzBhSU1NxdnamX79+LFy4kGLFit3155lbVE7GC+Unnp6e2qhrQ1q+hhYFmuznvE9rTcDcZ/CatJK61zSxtoqv+03m9YVTyclVsnKyj7XWLNrwBq37z6buNTPbqtlz+qNdDBpy59EZkTfJe7lwCAkJwd3d3egYhZqjoyNBQUF3LeoflDVv7X5TWFjYHUN+lFLBWmvP7NaVI9xCiAIlNTWeuS970vWNFdS9pgkp70DwD9sYvzhnxXZOKaUY9uQXJCyawtmSNrQ5F0+98c1ZsujOr1qFEEIUblJwCyEKjKSkq3zd24nRc/fjkASrHetRNfg8bQe3slqfjXpMxhy4iBNlTHhdTKTeW01ZHpD5GfRCCFEYRUREWPXodn4gBbcQokBISLzI8+NbMObnU5g0zGvyFF3CjlKlVnGr912zyWBMG78norSJphcTKfWmF7t3X7N6v0IIIfIHKbiFEPlefMJZWn/SgSXljjGqi4nZ3iN4Ye86itlZ55q7mXnEcyjJv37BxeKK9ueuM2F6N06fkdtTCiGEkIJbCJHPJSReoNtUb4L035Bsx5XGK3l507fYGPDb7bHWLxPx9Yt06W/PtmY7aT7+fZKScj+HEEKIvEUKbiFEvpWcHM2ip93w/yoct1N2DGQlyz/uSTb3Z7CqpkO+peuA7gCcqzeDQRN/My6MEEKIPEEKbiFEvpSaGs+8oR6MXH+eSnHwbORYlkzvYmixfdMrXRbTrZQrPcLNTFjiw5IFh42OJIQwkMlkws3NzfKIiIggMDCQ0qVLW67D3aZNG9auXZvp+lOnTuWzzz574BzDhg1jxYoVD9xOTjk6OnL58uVc6y8roaGhrFu3ztAMcqdJIUS+o3Uq373alFFL/wNg9uODeHXtZ3mi2AawsSmK//MB/DvTBdcLSYRP78CZjueoUSOPBBRC5Cp7e3tCQ0NvmRYREUHr1q0tRXZoaCg+Pj7Y29tne/ObvCglJYUiRYwrK1NSUrKcFxoaSlBQEF26dLmn9h7m9sgRbiFEvrNsdjeGfHMIk4ZvXDsyNnCJIWO276ZkqfrY/TSFG7Yw4MR5Phs8kgJ2nzEhxEPk5ubG5MmT8fX1zXT+gQMH8Pb25rHHHmP+/PkAxMbG8sQTT+Dh4YGzszM///yzZfkff/wRFxcXXF1defbZZ+9ob9KkSQwbNoy9e/fSq1cvAH7++Wfs7e1JSkoiISGBOnXqADB//ny8vLxwdXWld+/exMXFAWlHzF9//XXat2/P22+/TVRUFJ06dcLd3Z2RI0eS1c0V/f39cXZ2pnHjxrz99tuW6SVLluSNN97Aw8ODJ554gkuXLgHw77//0rlzZ5o0aULr1q05evToHf1PnjyZvXv30rJlS9zd3WnZsiXHjh0jKSmJyZMnExAQgJubGwEBAVy5cgUfHx9cXFxo3ry55Vb1U6dOZcSIEXTq1IkhQ4bkfOflhNa6QD2aNGmijbJ161bD+ha5R/azsbb98ZY+W1JpDdqvTmN9I9b80Pt4mPv4l2H1tQZ9vrjS/5sd/NDaFQ9O3suFw/79+42OoG1sbLSrq6t2dXXVPj4+Wuu0/39du3a9ZbmQkBDt5OR0x/pTpkzRLi4uOi4uTl+6dEnXqFFDnz17VicnJ+vo6GittdaXLl3SdevW1WazWR86dEjXq1dPX7p0SWutdVRUlNZa66FDh+rly5fr8ePH6xEjRmiz2ayTk5O1o6Oj1lrrN954Q3t6euodO3bowMBA3b9/f6211pcvX7Zkeffdd/WcOXMs7XXt2lWnpKRorbV+6aWX9Pvvv6+11nrt2rUasGS46ezZs7pmzZr64sWLOjk5Wbdv316vXr1aa601oJcsWaK11vr999/XY8eO1Vpr7e3trcPDw7XWWu/Zs0e3b9/+jv6vX7+uo6OjdXJystZa640bN+pevXpprbVesGCBpS2ttR43bpyeOnWq1lrrzZs3a1dXV8vP2cPDQ8fFxWW6H48cOXLHNCBI56A+lSElQoh8499zP7Pgq2/4LlazvXI52mzaR/ESeXuYRqf//cmerY40P5lAuVldiOwbSdWqeTuzEAWVtYadZfftVWZDSjJvJ+uGevbsib29Pfb29rRv3569e/fStWtXJk6cyLZt27CxseHs2bNcuHCBLVu20KdPH8vNZsqVK2dpZ/r06TRr1ox58+YBUKRIER599FHCwsLYu3cvr7/+Otu2bSM1NZXWrVsDcOjQId577z2uXbtGbGwsTz75pKW9vn37Ykq/je+2bdtYtWoVAF27dqVs2bJ3bMe+ffto164dFStWBGDQoEFs27YNHx8fbGxs6NevHwCDBw+mV69exMbGsmvXLvr27WtpIzExMdP+o6OjGTp0KP/88w9KKZKTkzP9We7YsYOVK1cC4O3tTVRUFNHR0QD06NEDe3v7LPfD/cpjX8IKIUTmomOO0Mb3VRZ4Xqdr19oU/XE/NWrbGR0rW8XsK1N20SRibaFfxAU+Gf6O0ZGEEHlUSEgIDRo0yHSeuu3TglIKPz8/Ll26RHBwMKGhoVSuXJmEhAS01ncsf5OXlxfBwcFcuXLFMq1169asX78eW1tbOnTowI4dO9ixYwdt2rQB0oZu+Pr68vfffzNlyhQSEhIs65YoUeKuOW93tw8VmW2z2WymTJkyhIaGWh5hYWGZ9j9p0iTat2/PoUOH+PXXX2/JmV2Gm7lv356HxdCCWynVWSl1TCl1XCk1IZP5ryuljiilDiqlNiulahmRUwhhrNTUG/Sd1ZtzxSIgpgo+g7bQrFP++XVQv+1ENvs8Qng5CCn1A79viTU6khCFktbWeTwMBw8eZPr06YwdOzbT+T///DMJCQlERUURGBiIl5cX0dHRVKpUCVtbW7Zu3crJkycBeOKJJ1i2bBlRUVEAtxTXnTt3ZsKECXTt2pWYmBgA2rRpw5dffkmLFi2oWLEiUVFRHD16lEaNGgEQExND1apVSU5Oxs/PL8ttaNOmjWX++vXruXr16h3LNGvWjD///JPLly+TmpqKv78/bdu2BcBsNluuovLTTz/RqlUrSpUqRe3atVm+fDmQViwfOHAg0/6jo6OpXr06AAsXLrRMd3BwsGzr7TkDAwOpUKECpUqVynK7HgbDhpQopUzAXKAjcAbYp5T6RWt9JMNiIYCn1jpOKTUamAn0y/20QgijaK1Z/F5bFvkeZZCPidq1lzFygKPRse7ZE9/+Rr2POhJZ4jxH5k0lss1nGHhCvxAiD9i+fTvu7u7ExcVRqVIl5syZk+UVSpo2bUrXrl05deoUkyZNolq1agwaNIju3bvj6emJm5sbTk5OADRq1Ih3332Xtm3bYjKZcHd3v6UA7du3LzExMfTo0YN169bRrFkzLly4YDmi7eLiQqVKlSxHfW8OQ6lVqxbOzs63FK8ZTZkyhQEDBuDh4UHbtm155JFH7limatWqfPTRR7Rv3x6tNV26dKFnz55A2tHlw4cP06RJE0qXLk1AQAAAfn5+jB49mhkzZpCcnEz//v1xdXW9o+233nqLoUOH8sUXX+Dt7W2Z3r59ez7++GPc3Nx45513mDp1KsOHD8fFxYXixYuzaNGi7HbVA1P3cmj/oXasVAtgqtb6yfTX7wBorT/KYnl3wFdr/fjd2vX09NRBQUEPO26OBAYG0q5dO0P6FrlH9nPu2r7uNRr2+ZLy8TDH5UnGhfxu9SuSWGsf+29/noFbfoBUWyaVO8C0VzP/6ljkDnkvFw4hISG4u7sbHUPkQMmSJYmNvb9vAGNiYnBwcHjIiW4VFhZ2x5AfpVSw1tozu3WNHFJSHTid4fWZ9GlZeR5Yb9VEQog85dzZDdi/MJfy8fBHtSoM2PBbnrv83714puUcOsTV4et1ydTyfYpMvm0VQghRABn5hWZmo+ozPdyulBoMeAJts5g/AhgBULlyZQIDAx9SxHsTGxtrWN8i98h+zh1aX+PCW8PoH5nMyZJF+O/l+RQL2w5h2a/7oKy5j1+r34oOn/1HEfNJRj/3BQNe8bBKPyJ78l4uHEqVKpXlEAiRt0RGRt73vkpNTbX6fk5ISLjv3xlGFtxngJoZXtcAzt2+kFKqA/Au0FZrnXj7fACt9TxgHqQNKTHqK0L5erJwkP1sfVqb+W5MQ14MiibRBOtHLGL0291yrX9r7mPdtg2/LlhLjx1X6B48g8cei6J6dblMoBHkvVw4hISEWH2ogTBebgwpsbOzu+/hSUZ+ObsPeEwpVVspVRToD/yScYH0cdvfAj201hcNyCiEMMCGbW/Sa9ExAL5qOoiRnw00ONHDo5QNrnNnElMUup2+ypcvfWB0JCGEEFZmWMGttU4BxgF/kPYl8TKt9WGl1DSlVI/0xT4FSgLLlVKhSqlfsmhOCFFAXLr6F/1+96dvX/ihbj1GrF9stZtVGOUR5+fY9GRVAHrt+pCjR1MNTiSEEMKaDD39SGu9TmtdT2tdV2v9Qfq0yVrrX9Kfd9BaV9Zau6U/ety9RSFEfpaaGkefb58j2u48W0vXw3FhMKVKF7Bqm7QbLLT66gcu2StaXIhn/htTjI4khBDCivLx+f5CiIJm2edPUSTsCKTaMqrCIrxblTQ6ktVUrNWZnV3TLsxU5vx3HD9uzCVahRDWZzKZcHNzszwiIiIIDAykdOnSuLu7U79+fdq0acPatWutmiMwMJBdu3bd13rduuXeeTR3s3DhQs6du+OUvzxPbrsghMgTjhzwpeUH2+l3HYa364/vpuZGR7K6lp9+Q9PS/dlX8wJ/zlrFprm9jY4khLACe3t7QkNDb5kWERFB69atLUV2aGgoPj4+2NvbZ3nzmwcVGBhIyZIladmy5R3zUlJSKJJH7saVmpqKyWTKdN7ChQtp3Lgx1apVeyjt5RY5wi2EMFx8/Cn+GzCBWtc1wRUcmPbdfAz+3ZgrKjl2xbV52p3hNqd8RESEHOUWorByc3Nj8uTJ+Pr63jEvNjaW4cOH4+zsjIuLCytXrgRgw4YNtGjRAg8PD/r27Wu5aYyjoyNTpkzBw8MDZ2dnjh49SkREBN988w2zZs3Czc2N7du3M2zYMF5//XXat2/P22+/zd69e2nZsiXu7u60bNmSY8eO3TVzQkKCJZe7uztbt24F0orinj170rlzZ+rXr8/7779vWWfJkiU0bdoUNzc3Ro4cSWpq2jksJUuWZPLkyTRr1ozdu3czbdo0vLy8aNy4MSNGjEBrzYoVKwgKCmLQoEG4ubkRHx/P5s2bcXd3p3nz5jz33HMkJiZafgbTpk2jVatWltvCG0kKbiGEobTWLHnJm25hN4ixheMTf6NW3WJGx8o17/ecQdHEUjyeHMynE78yOo4Qwgri4+Mtw0mefvrpLJfz8PDg6NGjd0yfPn06pUuX5u+//+bgwYN4e3tz+fJlZsyYwaZNm9i/fz+enp588cUXlnUqVKjA/v37GT16NJ999hmOjo6MGjWK1157jdDQUFq3bg1AeHg4mzZt4vPPP8fJyYlt27YREhLCtGnTmDhx4l23a+7cuQD8/fff+Pv7M3ToUBISEgDYu3cvfn5+hIaGsnz5coKCgggLCyMgIICdO3cSGhqKyWTCz88PgBs3btC4cWP++usvWrVqxbhx49i3bx+HDh0iPj6etWvX0qdPHzw9PS3tKqUYNmwYAQEB7Nmzh5SUFL7++mtLPjs7O3bs2EH//v1zuKesJ298dyCEKLR2bnyDPn7/AvB1yxd567XWBifKXVUrdOLriFI853+dX2pM5+LFl6lUyehUQhRM6n3rnIStp9z926nMhpRk2o7OvJ1NmzaxdOlSy+uyZcuydu1ajhw5wuOPPw5AUlISLVq0sCzTq1cvAJo0acKqVauy7LNv376W4RbR0dEMHTqUf/75B6UUycnJd827Y8cOXnrpJQCcnJyoVasW4eHhAHTs2JHy5ctbsuzYsYMiRYoQHByMl5cXkPZBpFL6LzyTyUTv3v8/rG7r1q3MnDmTuLg4rly5QqNGjejevfst/R87dozatWtTr149YmJiGDp0KHPnzuXVV18FoF+/fnfNn5uk4BZCGCYm9iiJo/9H2QRYX60aI9d8a3SkXKeUwvvtiSQsG0OPM5eZMm0l7/vKWG4hCqOQkBAaNGhwx3StNeq266NqrenYsSP+/v6ZtlWsWNo3hSaTiZSUlCz7LFGihOX5pEmTaN++PatXryYiIiLbG0Nl9QEBuCOvUgqtNUOHDuWjjz66Y3k7OztL4Z+QkMCYMWMICgqiZs2aTJ061XLkPKf9w63bZjQpuIUQhtDazMRZfXjvQiJXiin4dDOlyxS8SwDmRC2XkfzWZCLd9l6j5u9vEx/fG3t7o1MJUfBkdyTaSAcPHmT69Ol89913d8zr1KkTvr6+fPnllwBcvXqV5s2bM3bsWI4fP86jjz5KXFwcZ86coV69eln24eDgwPXr17OcHx0dTfXqaVdPWrhwYbaZ27Rpg5+fH97e3oSHh3Pq1Cnq16/P/v3yhOwCAAAgAElEQVT72bhxI1euXMHe3p41a9bwww8/ULx4cXr27Mlrr71GpUqVuHLlCjExMdSqVeuWdm8W1xUqVCA2NpYVK1bQp08fyzbcvIW7k5MTERERHD9+nMqVK7N48WLatm2bbW4jyBhuIYQh1ge/g2/yMRqPgU+6fchTA52MjmQYpWyoPmUoAIMj/uW7r0IMTiSEyA3bt2+3XBZw7NixzJkzJ9MrlLz33ntcvXqVxo0b4+rqytatW6lYsSILFy5kwIABuLi40Lx580zHf2fUvXt3Vq9ebTlp8nZvvfUW77zzDo8//rjlZMa7GTNmDKmpqTg7O9OvXz8WLlxoObLeqlUrnn32Wdzc3Ojduzeenp40bNiQGTNm0KlTJ1xcXOjYsSORkZF3tFumTBlefPFFnJ2d8fHxsQxBARg2bBijRo3Czc0NrTULFiygb9++NG/eHBsbG0aNGpVtbiOo7A7H5zeenp46KCjIkL4DAwOz/fpF5H+ynx/c1esHqftxZ64Wi6TkseeInPc9JfPQJbeN2MepqXH82agi3sfimNmwGW/+vQcbOSRiVfJeLhxCQkJwd3c3OkahsnDhQoKCgjK94oq1xMTE4ODgYNU+wsLC7hjyo5QK1lp7Zreu/DoXQuQqszmFZc934ZVdkdhedmTliC/zVLFtFJOpODaveAMw/MRefl513uBEQgghHhYZwy2EyFUblg/j2TVnKZ4CNhVG0KmddY9I5Ccths9j8+d1CK6cwNKAeTzdZ7LRkYQQ4p4NGzaMYcOGGR0jT5Ej3EKIXHM5ag9l31hG8RQIqFWPNxe/Y3SkPKWYXVWCvmjL250gpLI/hw8XrCF/QghRWEnBLYTIFWZzEmuG+NDsbDJnS5io+e1WuRJHJl5sM4ViSaWh4lHe+36T0XGEEEI8BDKkRAiRK35bNIBn/7gAgH+3j3nzyWoGJ8qbypVpQZ+y1XEIjKbayZe5cSOMPHQpWSGEEPdBjnALIawu8vwWar7zC8VS4cc6zrz845tGR8rTJnj5MPc3eCvsKAvmHjA6jhBCiAckBbcQwqpSUxMY8fVoLpRM4YSDLY0WbqFoUaNT5W0NWr1LYL3iFEuFa36vUsCu3ipEobR69WqUUtleKzs7w4YNY8WKFXdd5sMPP7zldcuWLe+rr6lTp/LZZ59lu1zJ+7jUVJcuXbh27do9r3ft2jX+97//WV6fO3fOclOcvEwKbiGEVf2w5QXW2oTTeYCJb19cRZPWFYyOlOeZTMVJfT7tsq6DI7axa+edtzQWQuQv/v7+tGrViqVLl1q9r9sL7l27dlm9z5zSWmM2m1m3bh1lypS55/VvL7irVauW7QeQvEAKbiGE1Zw+8xuvbV4PQJV/JjDjk24GJ8o/Wo7+mv/KmHC8bubXGe8ZHUcI8QBiY2PZuXMn33///S0F982bL/Xp0wcnJycGDRrEzRsSTps2DS8vLxo3bsyIESO4/UaFmzdv5umnn7a83rhxI7169WLChAnEx8fj5ubGoEGDgFuPQM+cORNnZ2dcXV2ZMGECAPPnz8fLywtXV1d69+5NXFzcXbfnxIkTtGjRAi8vLyZNmnTLvE8//RQvLy9cXFyYMmUKABERETRo0IAxY8bg4eHB6dOncXR05PLly7z99tu3FNBTp07l888/JzY2lieeeAIPDw+cnZ35+eefAZgwYQL//vsvbm5ujB8/noiICBo3bgxAs2bNOHz4sKWtdu3aERwczI0bN3juuefw8vLC3d3d0lZukoJbCGEVqak32DpwCOv8r1DnWCN+f28KReQ07Rwr4dCQkA5VAGh+eDFXrxocSAhx39asWUPnzp2pV68e5cqVY//+/ZZ5ISEhfPnllxw5coT//vuPnTt3AjBu3Dj27dvHoUOHiI+PZ+3atbe06e3tTVhYGJcuXQJgwYIFDB8+nI8//hh7e3tCQ0Px8/O7ZZ3169ezZs0a/vrrLw4cOMBbb70FQK9evdi3bx8HDhygQYMGfP/993fdnldeeYXRo0ezb98+qlSpYpm+YcMG/vnnH/bu3UtoaCjBwcFs27YNgGPHjjFkyBBCQkKoVauWZZ3+/fsTEBBgeb1s2TL69u2LnZ0dq1evZv/+/WzdupU33ngDrTUff/wxdevWJTQ0lE8//fSWXP3792fZsmUAREZGcu7cOZo0acIHH3yAt7c3+/btY+vWrYwfP54bN27cdRsfNim4hRBWsXJWNwbtuMLjp+GVMm/h2tjW6Ej5jsvbr5GioNvZiyz6OtjoOEIUDEpl/Zg37/+Xmzfv7sveA39/f/r37w+kFYX+/v6WeU2bNqVGjRrY2Njg5uZGREQEAFu3bqVZs2Y4OzuzZcuWW47cpm2G4tlnn2XJkiVcu3aN3bt389RTT901x6ZNmxg+fDjFixcHoFy5cgAcOnSI1q1b4+zsjJ+f3x193W7nzp0MGDAAgGeffdYyfcOGDWzYsAF3d3c8PDw4evQo//zzDwC1atWiefPmd7Tl7u7OxYsXOXfuHAcOHKBs2bI88sgjaK2ZOHEiLi4udOjQgbNnz3LhwoW75nrmmWdYvnw58P+F+81cH3/8MW5ubrRr146EhAROnTp117YeNjneJIR46CL+XYb7h9swafhfg1aM/WqI0ZHypTruY5nXega7yl9j5/4VvEoToyMJIe5RVFQUW7Zs4dChQyilSE1NRSnFzJkzAShWrJhlWZPJREpKCgkJCYwZM4agoCBq1qzJ1KlTSUi481yO4cOH0717d+zs7Ojbty9FsvkaUWuNyuTDwrBhw1izZg2urq4sXLiQwMDAbLcrs3a01rzzzjuMHDnylukRERGUuMv1Tfv06cOKFSs4f/685YOJn58fly5dIjg4GFtbWxwdHTP9GWRUvXp1ypcvz8GDBwkICODbb7+15Fq5ciX169fPdrusRY5wCyEeqpSU6+waNJLHrpo5VNaeDsv+wGQyOlX+ZDLZYTPJmyWucKJGAKEHzEZHEiL/0zrrx4gR/7/ciBF3XzaHVqxYwZAhQzh58iQRERGcPn2a2rVrs2PHjizXuVlYVqhQgdjY2CxPCqxWrRrVqlVjxowZt9xK3dbWluTk5DuW79SpEz/88INljPaVK1cAiImJoWrVqiQnJ98xDCUzjz/+uGUsesbln3zySX744QdiY2MBOHv2LBcvXsy2vf79+7N06VJWrFhhueJIdHQ0lSpVwtbWlq1bt3Ly5EkAHBwciImJuWtbM2fOJDo6GmdnZ0uur776yjIOPiQkJNtMD5sU3EKIhyrgwycZ+Nc1km1g+3OLqde4uNGR8rV+Xu9gn1QWyp5g+uLNRscRQtwjf3//W05uBOjduzc//fRTluuUKVOGF198EWdnZ3x8fPDy8spy2UGDBlGzZk0aNmxomTZixAhcXFwsJ03e1LlzZ3r06IGnpydubm6WS/5Nnz6dZs2a0bFjR5ycnLLdptmzZzN37ly8vLyIjo62TO/UqRMDBw6kRYsWODs706dPn7sWxzc1atSImJgYqlevTtWqVS3bFRQUhKenJ35+fpZc5cuX5/HHH6dx48aMHz/+jrb69OnD0qVLeeaZZyzTJk2aRHJyMi4uLjRu3PiOEz1zg7r9rNf8ztPTUwcFBRnS982zjUXBJvs5a+Fh31O0xYs4Rmtmu3TkpZAN2OTDj/V5bR+PnFwfj3XhJMU8xsi/w+U65g9JXtvPwjpCQkJwd3c3OobVjBs3Dnd3d55//nmjoxgqJiYGBwcHq/YRFhZGgwYNbpmmlArWWntmt24+/FMohMiLkpKi6OX/IV97aXZVKkXPlb/my2I7LxrZvBMjg2H4iX9Y7p+7J/oIIfKuJk2acPDgQQYPHmx0FJEN+XMohHgopq8eyGHTf8z0Ksn2qftwfLRY9iuJHHHuOJW/ahajZDKEzXvb6DhCiDzi5mX3Mp54KfImKbiFEA/scPAsFqZfO7bR2Q95a1Q9gxMVLLa25TnbzRGAJ06sIzLS2DxCCCHujRTcQogHkph4gRNDJnFg/g2e3OXBho/G3eslakUOeL06kfgi0D7yOot9txsdR4h8paCdryZy34P+H5KCWwhx37TWLHu9Pd2O3KBoKgzsOJtq1aTatobqjw5gW/2069gmr51qbBgh8pHU1FSioqKk6Bb3TWtNVFQUdnZ2992G3PhGCHHfDuycTpcFYQDM8RrEO2+0MjhRwWVjY4vu1xgm/4XPmR0cOWymYSM5ZiJEdm7cuEFMTIzlFuiiYEpISHiggjg7dnZ21KhR477Xl4JbCHFf4uJOEjX8Y9ziYUO1iryw6kcZSmJlXiM+4pPfO+PvnIT7T3+x4IMWRkcSIs/TWlO7dm2jYwgrCwwMzNOXf5TDI0KIe6a1mYDR7XnieDxXiymuTPqdSpXl14m1lavUjg09a3OgKqz4ZzFmufGkEELkC/IXUghxzzbtfIuey04A4Pv4aPqP8jA4UeGglOLldm0BiHVcwbYdd966WQghRN4jBbcQ4p5EXQvmmd9+ossg+Kp+A15Z5Wt0pEKlQ6M3GBpaih1LL7Fh1odGxxFCCJEDUnALIXIsNTWBAd8N45pdJH+VqIPTvL8oVVoGbuemEiXq0S3JnsdPQ8PQBSQlGZ1ICCFEdqTgFkLk2C9ze1Is+BCYTQwttYCObRyMjlQo1X3pGQB8zp7k55UXDE4jhBAiO1JwCyFyJCL8J1ynbOLXpfDclp7Mn9zG6EiFViPv9wiuVpSSybB//jSj4wghhMiGFNxCiGwlJUUR+swo6lwzE1K+OG/+bxG2tkanKryKFq3EiTaVAfA8vpobNwwOJIQQ4q6k4BZC3JXWmuVvtMXnQAxxRWDPmGU0cC5pdKxCr8HY4QB0iYxk5dLTBqcRQghxN1JwCyHuateG8XSdfxiALzz7Mer9rgYnEgD1mr3KXzWKYp8CRxdPNzqOEEKIu5CCWwiRpcuX92DzwhzKJMJv1asy6pef5G6SeYStbVkO9HVkWE/4pvIJYmKMTiSEECIrUnALITKVmnqDsTMHUzUmmbMlTNj77qBCRfmVkZe0fOklFrnZcLV+IP5rrhgdRwghRBYM/euplOqslDqmlDqulJqQyfw2Sqn9SqkUpVQfIzIKUVh9uLofy0r8i8eLRVg45Hu8feoYHUncpn6NwdRMqQGmFP63ZZnRcYQQQmShiFEdK6VMwFygI3AG2KeU+kVrfSTDYqeAYcCbuZ9QiMLrr4MzmRoaCLZQ4fwHvLNoqNGRRCZsbcvwfJUKVFp0iuJXphEdPYrSpY1OJYQQ4nZGHuFuChzXWv+ntU4ClgI9My6gtY7QWh8EzEYEFKIwuha1H7pO4rPNNyhx/Cm2fToeGxlJkmcNbDmAEcEwMCKS5YvDjY4jhBAiE0b+Ga0OZLyW1Zn0aUIIg6SmxrG+z5M0O5PEM3+bWNp1FlWqyFmSeVltlxHsdrTD1gyn/N83Oo4QQohMGDakBMjsr7i+r4aUGgGMAKhcuTKBgYEPEOv+xcbGGta3yD0Fdz9rjq4aw6jAy6Qq+KL9e3SvEklgYKTRwXJdftvHZ1tXpNWJ0zwesZ5ff92Bg0OK0ZHyhfy2n8X9kf1cOOT1/WxkwX0GqJnhdQ3g3P00pLWeB8wD8PT01O3atXvgcPcjMDAQo/oWuaeg7ufA1aPpN/8oALNduzBzzVRMJoNDGSS/7ePjZV4jZfHreJ+/yqKIqnR/qa7RkfKF/Lafxf2R/Vw45PX9bOSQkn3AY0qp2kqpokB/4BcD8whRaJ36bxUVR35H2QT4pUYNhm74pdAW2/mRY+MX+KtWMWzNELFshtFxhBBC3MawgltrnQKMA/4AwoBlWuvDSqlpSqkeAEopL6XUGaAv8K1S6rBReYUoqOITzrB5wHM0upRCWJliVP5xH+UrSrWdnxQp4kBk2yoAeEWsJzbW4EBCCCFuYeSQErTW64B1t02bnOH5PtKGmgghrMBsTmLwtz3Z0j4a++u2xA5cwwvtqxgdS9yHhiOe49PwKSx3iuLa2usM7V/K6EhCCCHSycW+hCiktNbMWNmHVdf2c81O8WOnpbwwqbPRscR9esxrLLPaV2dfrRS+2Syj84QQIi+RgluIQmrj0uco8+k6bFPA8fg01nzay+hI4gHY2panW920byf23VhNQoLBgYQQQlgYOqRECGGMI0Ff8eiYH+l0zcz1RGdGbX6XokWNTiUe1Aste1L09WC6hK/h99+u4tO7rNGRhBBCIEe4hSh0LkVu5/rTb1LnmpmgCiXo+V0gFSrIzW0KAtfaz/NiSBG6/Gtm3/efGB1HCCFEOim4hShEEhLOsffJrjQ/k8TpkibOfrYTZ69yRscSD0mxYtU45lUBgHphy0mR+98IIUSeIAW3EIVESkoMAd096Pp3DDFFYdXIJfQc6mp0LPGQ1Xq+KwDdz59gy6YbBqcRQggBUnALUSiYzUn4jvVi6KYLpCqY9dQkXvmsv9GxhBW4dJrI0fImyiVots3/0ug4QgghkIJbiAJPazMT/XsyscIxfnsMPm8+kPdWTTM6lrASe/s6HGySNkyo1oElmM0GBxJCCCEFtxAF3dz1w/jk2Gbii8Irzd7nla1+2Mg7v0CrNrQdAN0jw9m9O9nYMEIIIaTgFqIg++O7QTi8408xczIV/hlDsO8kihUzOpWwNrfuk/jRpRiTnjDz9ertRscRQohCT67DLUQBtSlgNO6v+FMpTvOfqRGj1n1F6dJy+b/CoETJxiwa5MiW+GOUOfwrWnujZNcLIYRh5Ai3EAXQn6tepuEL31IpTrOxWkUG+f1F1Srydi8slFIMbe4OwLXK6zh2zOBAQghRyMlfYCEKmD2/vsWjQ32pFqsJrFKW6muPUq9BCaNjiVzWpfEI2v5nx+y94aya/7vRcYQQolCTgluIAmTHz29TY9BnVI/VbKtSmnI/h9PQXW5sUxiVL9uG8aFFeHkv6C2zjI4jhBCFmhTcQhQQP++ZypH3v6RGjGZ7lVLYrziGS9MKRscSBlHKRErn2gC0jdzNhQsGBxJCiEJMCm4hCgD/bW/S65cveLlLEp8716Lk6n/weryy0bGEwVyff4VEE7S8EMOqxYeMjiOEEIWWFNxC5GNaa1Z+1ZsX1n6NuVgM6mxXnlpxDPfmlYyOJvKAmnX7s9vRDhsgcvUnRscRQohCSwpuIfIpszmFH8e3pMvrq1i+Oo4yx304PHU1DevJhbZFGpOpBJdaVwHA89RGbtwwOJAQQhRSUnALkQ8lJ8fwYz8nBn+xB/sUuIwThz4MoE4tW6OjiTzG6YVBAHQ4f4H1v1wyOI0QQhROUnALkc/ExZ5iddu6DFvxLyYNvq5teTr4MNVrFjU6msiD6nmO4xcnO+Y1gSUbfzU6jhBCFEpScAuRj5z67zd2uDfkmd2XSLKBDzo8z6igQBxKyVtZZK5YsSr4v1SP156C3+MCSU01OpEQQhQ+8ldaiHzi9+BP+WFYfzodv8EVO8WXQz9h4obvKFLE6GQirxveoj0AiY+sZ8cuqbiFECK3yZ9qIfI4szmFGSuGMvXgGoq0jeORS2UpOTqAt17uaHQ0kU+0fOxF6kT9SIcTl1l7+Ufath5udCQhhChUsj3CrZQqrpSapJSan/76MaVUN+tHE0JEXz3MN51r82XwT2jbOIqdeRqPFWd4RoptcQ9KlGjIJyHw7VqoucvX6DhCCFHo5GRIyQIgEWiR/voMMMNqiYQQAOxZN43wRh6M2XiGRasVDU9P4/Sslbg1Km50NJHPKKUo1ccNgCcv/E1YmDY4kRBCFC45Kbjraq1nAskAWut4QFk1lRCFWHLSdRY/14RGT0/BKzKJUyVNnO70JYfmT6JMGXnrifvjMeAdouyh/rVk1sxba3QcIYQoVHJScCcppewBDaCUqkvaEW8hxEN2aM98NjjX4NkF+3FIgpWO1Tj783HG+L6MklpbPIByFdqzu17JtBeBc4wNI4QQhUxOCu4pwO9ATaWUH7AZeMuqqYQoZJKTY3jrq55U8R5J1/AYYorCp50G8+TB07TwdjQ6nigAbGyKkPJUHQDanPuLCxcMDiSEEIVItgW31noj0AsYBvgDnlrrQOvGEqJw0Fqzdt+XVHu3CZ9e+YXFrpo/alRg+3e7Gf/HYko6yJU7xcPjNvwlEk3Q4mIMa346bHQcIYQoNLK8LKBSyuO2SZHp/z6ilHpEa73ferGEKPjO/beVwMHP4v9IJJcbmCG2En80e48RG8dSoqQU2uLhq1GnP1vrvIxOjeeP7csZ+VojoyMJIUShcLfrcH+e/q8d4AkcIO1kSRfgL6CVdaMJUTDFXP+PgFF96P5zKAPjNF7HIKjYEJaN+4LWXuWNjicKsCJFSvLbu42ZHbEPm/D9xMeDvb3RqYQQouDL8jCa1rq91ro9cBLw0Fp7aq2bAO7A8dwKKERBkZR4hZ/GtifSsR4v+IdQOU6zq3JJtoycz7mFi6TYFrliUAsfAMy1N7F+U7zBaYQQonDIyffWTlrrv2++0FofAtysF0mIgiU+8SLTP+3H0ZpVGPi/QOpdTeV4aVtmPz0S5/BoRn74glyBROQa50eGUDG+Ci0j4/l96fdGxxFCiEIhJ7d2D1NKfQcsIe3SgIOBMKumEqIAiL7+L68tnozfya1QLJLnEuBMSRuWNu/J4IV+vFJdvssXuc/OrgZfHtYMXAvf1ZmL2TwOGzllQAghrConv2aHA4eBV4BXgSPp04QQt9Fas3/rfL5v9xiXa9VjzamfSCoRSXJMXaY//TKpIVd4c+MqqkixLQxUo3/aKThPXQwnaF+KwWmEEKLgy/YIt9Y6AZiV/hBCZCI+7gK/TB9HmeUb6fBfNB7pd87uv6ca59u+y7xJI6lQ3mRsSCHSuXWbwBmHVdSIMfP9PD+aNhtqdCQhhCjQsi24lVInSL/LZEZa6zpWSSREPpGSEk/Atm9Jenc2HQ+dpF9s2tskyQZ+cazGlR7j+OKjt7Gzl+/rRd7iUKoJmxqUosbeaErsmQdIwS2EENaUkzHcnhme2wF9gXLWiSNE3paYcIXwjd8w8Y+Z7E/5m8SSZ/j9EtSIheOlbfmtcXPc3vmKp7u6Gh1ViCwppbDzaQB79/DEhf2cPAm1ahmdSgghCq6c3GkyKsPjrNb6S8A7F7IJYTitNSfCtrB4dHdWNa7CtfIVGPFhALbn15NY8gzcqMQ3Lk/x43uzqRGZwCs7ttFWim2RD7gNfpOYouAWlcDPC3YaHUcIIQq0nAwpyXjHSRvSjng7WC2REAbSWhNxYT/f/7aS+t8updHJc7heSuTZDIOqIhxMNIxsTfOmY3ir/9OUL5uTL4qEyFuqVO/GxrrFaXwmjqBdfsDjRkcSQogCKyeVwucZnqcAJ4BnrBNHiNx1+fRB/lzyA1d27MYm8iyjvRXJDmcAOBEOjtGQbAM7q5Qi9FEXqvYdTZn6Vfm6U3uDkwvxYGxsirFvvAudT+5BnfgH3+tQqpTRqYQQomDKScH9vNb6v4wTlFK1rZRHiIdOa82Fq//xe9B2wv/8k0e37qbSxUs8duU6j11NoXeGZcd3hKvJdthdcWbO46VpVL8pHUe/SpvHKtImfZnAwEADtkKIh69X+8FMWrgXXetPVq+LZmj/0kZHEkKIAiknBfcKwCOTaU0etHOlVGdgNmACvtNaf3zb/GLAj+l9RQH9tNYRD9qvKFhSUpM4Gb6f//7azfmDB4n/7wRFz0dS+uo1KsfEEtBQMaf1DQA6noIPd///ukk2cLi8HWEVKnO1dkM+aTQCn65dqFi+qEFbI0TuebTaM1RJ+pjklDP8+stChvZ/xehIQghRIGVZcCulnIBGQGmlVK8Ms0qRdrWSB6KUMgFzgY7AGWCfUuoXrfWRDIs9D1zVWj+qlOoPfAL0e9C+Rd5kTknm8tVITkZd4PSlC1z5919Mhw6TeOUK5uhr2F67hn30dYrfiKVUXBw9fKqQWOwaZrsoguen0PF85u2eqACk2lIkug5nTFWY3/w6qbWdqNasLV69++FeowzuubqlQuQNRYtWZMoFG15cAEsfmUVKyisUkVMShBDiobvbr9b6QDegDNA9w/QY4MWH0HdT4PjN4SpKqaVAT9LuZHlTT2Bq+vMVgK9SSmmt77guuJESk5IICTvI2UMh7E9JwIwZ0iMqnTakIbFsWcx2aZ9TikRfw3TtetrKWoPWaDRo0CYTCTWqYdZmtNbYnYhAp6SSajZjNqeCNpNq1mizmYQypYkvU5pUsxlTzHVKnL+AOdWMTl/XbDajzWn/XnZ8hBSTDWZtpvSZcxSNjbXM01qjtRnMmhsl7LlcrQpms8aUmECtQ2GkJiVhTk5GpyRjTk7BnJqETkkluE4NzjsUJzk1FaeI09Q/HYlKNWOTmkyR5BRsU1IwJSdz3VYxvWVtknUyySSx9NcjlEtIplhqKnYpZuxTzTgkmSmVCB93gln/1959h0dV7V0c//7SCCEkhBaaFJFmoWhEEZEAgihVLBRR9KrYUNGLCpYLgly7YveiKGABQZAiKAiKqKgIgiAgRWpooUMCaTP7/SPjlZcLGEwmJ2V9nifPnJnZzlmTbWBxcubsi7K/NT2Ww7hJJ/6+R0Ye4EhU9vbGmHBKpsO26CiSS5fhUFwFfJWrU7puA6pfmMjBVomUjlaTEDlWk65tCX1nFJfv3Mw3X6fRqk2uj6eISBHh/P7s/uEcPr8Pf1YWvowMfFlZ+P0+fFk+fL6s7HE+H74SEfjN8PkcHE6FzAz8vuzu4vyBLwcuIhx/6Wj8zoHPR9j+/dldhEAnCrQ85xyZpUrjSpQAHCGphwlJzf6NtSO7vwD4S0SyY/9uj75LOXPCBuKcmwpMNbNmzrnvTzQuF6oCW466nwRccKIxzrksMzsAlAMK1Hf1x5W/0HtMUzaPOPGYy6+Dz+tkbw+bC49+c/xxW2Kg+v1/3k9+BiocPv7YR1vB8JbZ211XwScfnXj/FQfArjzQO4MAACAASURBVOjs7Znvw+Xrjj9ucn24s0f2dvX9sOnlE79m++tgVh3AoMV6uOME72lzDNzR6c+PATTZBRVP8J6i0sOw1HKEZkazE+OrartIjYjgSIlIDkWV4khMLP64ikRUrsrwJi2pU7s+DU+vS4V/RWEGDU4cV0SO4+zWA1hTdjR19/qY9/abtGrT3+tIIkWa3+fj0L7tbEvZz7a9e9ixbze+pSvwb00ic/9+fIdSsCMpWHoGlpHB5tiSTKtXjXRfBlGHUxn09QoisrII9/mI8GUR4fMR5neE+f0MahHPd1Ujcfi5/Ze93PnLAcL8jlDnArcQ5nfsijLO6lsKzAchWWx8JYMqhyDEZR8o/OOa0Qb8+xIYHLgYdIfV8Om4E7+3qvfDtsCHr6d+CJ3XHH/cjDrQ8brs7fhDsOP5448D6NwDptfP3n5sHgyd979jRjWBIc3a0aPr1Sd+IY+d7JSSB51zzwC9zKznsc875+7J5b7tOI8de+Q6J2Mws75AX4D4+Ph8/1Db6m0b8aWVJyl6b3bAQOo/gjogIzOWkJQSAOwPSWV97OE/n7c/x24vFUrIofKBe8aaMrvYHekPPG+4o15/r8USvj8WnJGaeYRfyu/6nzEOwxmEHahBZHoYYPxeajs/xR/GBb69f+YwNkRFE72rWvbejmTxWY1N+EMMn4XgtxB8ISH4Qgy/heLPrEnNXWUIJZT1ZXbzTpNd+ENC8YeEkBkeTmZYGL6ICNIjo+iUcRYlwyMoGR7B2D7JlIiIIDwqirBS0ZQoFUuJMnGULFuBNqWjuPToWe8H0WR/VTjeNz91PyuXLzyV6cq1lJQUfXCyiCtec+zYeXYMdefvo/zi0cyb19jrQPmmeM1z8RXseXb+TPZvXcnutWtZm5HBynA/e9NSid+2nR5L1hBzJI3Y9HTi0jIpn5ZFXJoj1kGD+2F7oJxO+xA6naCczjwDngh8aq5iKFzz24mzRIdtJD0uezvGoN6+E2Q2ByUO/fd+uA/C/Sd7kwbO8AOZIf7/9oxjb8ksCRmh4CAlLI39JbICz9n/G3swNBw7nH2FaZfmSI76/0GPLnqZWTFYavbnqg5zhB1Rfx6x+6PrHAyJJNIiC/TPs53o7Awz6+Scm25mx13z1zk3Jlc7NmsGDHHOXRa4Pyjwuk8eNWZWYMz3ZhYG7AAqnOyUkoSEBLdo0aLcRPvb5s2bR2Jioif7lvyjeS76itscz36tLe36zWFdbBhZ36dTv8FfrolWJBS3eS6ucjvPvqwMViUtYe6Sn1m0YS3Nps6h6o6dVDyUQqXUdKqk+ijhyx775MXw8KXZ263Ww5djj/+aaaFw9j9i2VA6hlBfSQZ/f4Bm29M4HB7OkYgSpIdHkBkegS88nK3ly/JVowaUjChBaUJovWIdISVLEhIZSWhUKUIjIwmNiCQ0IoKU6qfh4uIoER5O6ZTDRKelER4RQXhkScJKRFAiMoqIyJJElIwkrHx5QkPCiAgPJyzLT0hIKGHhYYSEhBIaFkpoaBghoaGEhBh2vMOfBYxXP89mttg5l/BX4052Ssn0wG2uivVJ/ATUCVxicCvQA+h1zJhpQB/ge+Bq4MuCdv62iEhhd961D7H7gTmccSCLZ96eSf3nO3odSSTf+X2ZrF74KT9PmcaRJcsptyWJ2rv3U+lwOuc8wH9/537/GmhyzIf095UwtkWHc8RXmaq7ziY2PJbo2DBe75JMWIV4SlWuSlzNWlSuXY9q9c6iXMWyrDuFf9cOzrN3KV452Skl0znO6Rt/cM51zs2OA+dk9wNmkX1ZwHeccyvMbCiwyDk3DRgFvGdm64C9ZJdyERHJQ3HlEvm0bmk6LDvE9gXvkv15eZGiyznH5h1LGPPVDGavWs7pSxbz4twNNEhzx/0sUOXkMuwNqUGcq8JH56WxoGQsUafXo/LZjWnQrAXVa1clzmAo2V8ixzrZZRueC/bOnXMzgZnHPPavo7bTgGuCnUNEpDgLCQlja79ziF+zgD17d/DwLqhw3A9NiBROzjk2rpjN/JdfpeQPS2ictJOPzsnK/jBgGKRXgXJpsDfS+LVsNBvKV+JQjTrEJjSnYbsubDn/LEJDvX4XUpid7JSSr//YNrMIoD7ZR7xXO+cy8iGbiIjkk8va38Zdm3+Cat8zfnoyd/+joteRRHIlPWMvM18dzMGPPuWXjdtpnJzO0ctkn59klEg+j1qhZ9LkjEbMmdiQFh1bc0lkyH9XFhbJK395YWIz6wC8CfxO9hlMtczsNufcZ8EOJyIi+aN6pa5U9z/GZtvMx1+O4+5/aNVJKXwO7d3Aa7PeYszS71ljyxk7cw99lmc/lxYKCyrFsqrWWZRqfSXtbr2dtGrR3gaWYiMnK4E8D7Ryzq0DMLPawAxAhVtEpIgIC4uhv68k17wMC2OfJi3tXiK1Bo4UAlmZh5n12gMcfncKbVdvY2EX+O3M7OfGnhFPKnFkXdyJ9v0foHXdCrT2Nq4UUzkp3Ml/lO2A9UBykPKIiIhHEi/vRNWnVtMubTtzPt9Px65lvI4kckK//TiOnx4dxkWL1tBhv++/jzddV5518d3p26IHfQc1Z8GCr3X5R/FcTgr3CjObCUwg+xzua4CfzKwbgHNuchDziYhIPqnftB9LK42g8Y4sFr77Ih27Pu51JJH/x+/P5MP5L5L2wAj6/Lyd+oHFWpKiQ5hR5yxirn+IAXdfx8CctBuRfJST/yUjgZ1AYBFxdgFlgU5kF3AVbhGRIiAysgarGpWl8Y5kqv86AeceLxQLXkjRl562m+HjH+PlX+dwoPQ6BsZmL0M+tWZFNlx6I9cNH8ptFUt4HVPkhP6ycDvnbsqPICIi4r3K1zWHWZ9w+c41/Lw4i/MSdKhQvJN+OJnJ911L448XkH5uJgcuBo7E8flpnWg2ug+dr2+tfxRKoZCTq5TUAu4Gah49PrcL34iISMHTpNNAkkpPodohP6NHfsB5CX28jiTFUGb6PibdeyUJH31Hz/1ZAFy1PJwF5zzMf+4YwJl1dHURKVxycuhiCtkrPk4H/MGNIyIiXoqJTWDOmTFU+/EAYd+PAlS4Jf/4/VnMePYGTntxEj12Zi/5sS42jMkXdKHHyNF8U0NFWwqnnBTuNOfcy0FPIiIinjMLIePWJlx01jx+DD3IdUlQrZrXqaQ4+HbFe7z2xBDGjV8PZH8Q8oOLutFr5FgerFHS23AiuZSTwv2SmQ0GZgPpfzzonPs5aKlERMQzLTvdy/Ubf8Af9gtjpmzmkX7VvY4kRdiBQ+vo/cYdfHroG6iXzs01Q1he4SLajZzEQ4214qkUDSE5GHMOcCvwFNmL4DwPPBfMUCIi4p34cpdR2yoD8OHCiR6nkaLKOR8znurB5toNWLF1DoSlE7u5Gwff/J37Fn7DWSrbUoTk5Aj3lcDpzrmMYIcRERHvhYaW5LYyZTjnNUhPf4qUlH8SrVNnJQ/tSlrAN12vpNvi7HX0BnwVy7LrxvLaW50JDfU4nEgQ5OQI9y+AlhsTESlGLr/kKlpvgHbbdzNjUpLXcaSIcM4x66XrOXT2JXRbnExGCDx/bnPaTtzCmw+pbEvRlZPCHQ/8ZmazzGxa4GtqsIOJiIh3aje8lYXVIojww8oPdBah5F7qwfVMaF2dtv3f5/QDPpaUj+SdR8Zw/6JvqVOvtNfxRIIqJ4V7MNmnlfwbeAFYCJwRzFAiIuKtiIiKbDy/PAB110zF5/M4kBRqyzZMpOO9rbhyfhJ+gxcbnk+pb5K5fegNWrhGioW/LNzOua+BA0AHYDTQBngzuLFERMRrNfpcBsDlOzfx7fw0j9NIYeScn1dm3sm5I/syr+Zmbm8bw7M3P8k9Py+kbn0d1Zbi44QfmjSzukAPoCewB/gIMOdcq3zKJiIiHmrUZgBryo6l7l4fX7/1Ji1b9fc6khQiGWm7mdw1gVllN+GrBxFJLenywAd0aVPV62gi+e5kR7h/I/todifn3MXOuVcA/VJRRKSYiIpqwJJzYgEot3isx2mkMNmz7Ue+Orc2PWZt4r1P4Ozf7mbjsLkq21JsnaxwXwXsAL4ys7fMrA2gM61ERIoJMyOqTwK3doInLt7HmjXO60hSCKxe8AZbm7TgslUHOVACXu4wgKXvv0zlSroEiRRfJyzczrlPnHPdgfrAPOA+IN7M3jCzdvmUT0REPHRh5wcY3TCaHdU3Mmrqaq/jSAH3zfv9iGl3Fw2TM1lbJoz3HhzP4I+f1eX+pNjLyYcmU51zHzjnOgLVgKXAwKAnExERz5WLa8mZ4fEATPxlksdppKByzjH9yY40uvk1Kqc6vqwSw9q3l9FvaHevo4kUCDlZafK/nHN7gf8EvkREpIgLCQnnplqV4fXfqbjrJfbseYRy5bxOJQWJcz6GTOjB3N9mcakfJtSqytlTV3LmOTFeRxMpMHJyHW4RESnGulxwHXcsgh7rdzH1w1Vex5ECxO9P545R7Rm6ahLfnZ5Fu/ZXcvE3G1S2RY6hwi0iIid1Wu2eLKgZSQiwZeIzXseRAsLnO8JH3Rqw9cs5YI4a6x7l0/cnUaVquNfRRAocFW4RETmpsLBYtl9UEYCGG2aRkeFxIPGcLyuVSZ3q0XPqBiZMhMTfHuW3t4YRG6uLmYkcjwq3iIj8pQb/uAqAdju3M3fWfo/TiJd8Wal8ckVdrv1sCz6Dfzfvzeyxw4iM9DqZSMGlwi0iIn+p/gV3s7RSGKUy4cd3RngdRzzi86UyuUM9rv5iG1kGjyXezOA57xGus0hETkqFW0RE/lLJkrVY2agsAKct+winNXCKHb8/nYlXn801s7dml+3WtzN09tuEndL1zkSKJxVuERHJkfgbmjPiAnjn/N388osad3Hi92cxeEQbOs3YCMDjl/Rm2OdvqGyL5JAKt4iI5EhChwcZ2KYMCxrs5j/TfvY6juQT5xyDJ3TjiQPfc1lveKRZNx6b/Z7KtsgpUOEWEZEciYlpSuOSFQCY+ttkj9NIfnDO8fqEa3li1WcQ4mdT5iAemTOJiAivk4kULircIiKSI2Yh3NCwNrcshkHLRrJtm9eJJNjmvtGDnjd+TLffsii/5m6WvTycqCivU4kUPircIiKSY50b9eGFWXD3it1Meed7r+NIEC357FEaD5hI2TS4Ykltfn3hJeLidJ1tkb9DhVtERHKscrXOfFs7+xDnvunPe5xGgmXT8g+J6fUU5Y84ZlYrS/Opy4iPV9kW+btUuEVEJMdCQ6M4kFgFgKZb5pGa6nEgyXP7dvzEnvY3UXu/j0UVIokatYz6DXQeiUhuqHCLiMgpafiPPvgMWu7cw8yPk7yOI3ko/cgOFrduzbnbMtgQE8qKx+eT2K6q17FECj0VbhEROSVnNLyVhdUiiPDDirFPex1H8ojfn86gxy+l2boU9kXC+Js/pM8d53sdS6RIUOEWEZFTEhERz4Zm2ZcHPHvtFNLTPQ4kueac4+kp3XkxciUX/QOevPxfDHz+Wq9jiRQZKtwiInLK6t12LV/WhJkNk5kx+4jXcSSXZv34KI/+/AWY40Dqvxj+0eOYPiMpkmdUuEVE5JSd3XwA3XtW5t3zM3hl5gyv40gubFo2ltpXPE3vFYcp+ftVLHx+COHhXqcSKVpUuEVE5JSVKFGFy6pXBOC7fZPJyPA4kPwtKft/I7lDX+rs83Hb91HMufNdKlbUoW2RvOZJ4Tazsmb2hZmtDdzGnWDc52a238w+ze+MIiJycre16MgFW+CJjR8zZ8Zer+PIKfL5DjOnUwvOT0onKTqElYPmcdGFpb2OJVIkeXWEeyAw1zlXB5gbuH88zwLX51sqERHJsfNr9+X1meE8+GMmP731jNdx5BRNe6glXb/dTUYIvNn1eW65V1ckEQkWrwp3F2BMYHsM0PV4g5xzc4FD+RVKRERyLjKyOqualgOg3oqPyMryOJDk2M8zHuTSVxYB8O9zL+fx0f09TiRStHlVuOOdc9sBArcVPcohIiK5cHrfTgBcsWMjX80+6HEayYnd+37C13cEpTNgYs3K9J06ndBQr1OJFG3mnAvOC5vNASod56lHgDHOuTJHjd3nnDvRedyJwADnXMeT7Ksv0BcgPj7+vPHjx+cm+t+WkpJCdHS0J/uW/KN5Lvo0x6diG2V69qHxjiweaHkvHYYc9xeWBVLxnOfDDPzmPg7uWMPwWVH8fvPHJDQr6XWooCqe81z8eDXPrVq1WuycS/ircWHBCuCcu/REz5nZTjOr7JzbbmaVgeRc7mskMBIgISHBJSYm5ubl/rZ58+bh1b4l/2ieiz7N8an5IOF+Gn+6k8Ybp9OixYhCc7S0uM2zc45B49rxo38NlC3BnBu/5LVBF3gdK+iK2zwXVwV9nr06pWQa0Cew3QeY6lEOERHJpZq3tAeg4/b1zJ+b4nEaOZGFE+9i5/h54ODMTc/yysNFv2yLFBReFe6ngLZmthZoG7iPmSWY2dt/DDKzb4CJQBszSzKzyzxJKyIiJ9Sk7cPMqRXB+43g7claBKcg2r11PrF3vsW707O4/fPz+fqFfoRoJQ6RfBO0U0pOxjm3B2hznMcXAbccdb9FfuYSEZFTFxVVl9furcuU/b8SufJzsrK6E+bJ3y5yPD7fYb7v1pVOe7JYGRfONcOmUb68FrcRyU/6962IiOTaXYltAUirNYXZX6Z7nEaONvPR1nRauI/0UJjUYySt2x7vegYiEkwq3CIikmsX17uXGgcqcNOq/cx+/Vmv40jA6u+e5+IRPwLw7yaX8/ArN3obSKSYUuEWEZFci4yswfDN4bwzDS5Z/CbpOsjtuSMpG9jf8xHi0uDT08px00Rdb1vEKyrcIiKSJ+r364kfuGL7Vj6dsNXrOMWac34Gvnol0WnpbCtl7B40l5o11bZFvKLCLSIieeKcix7gh+oRRPpgxdtDvI5TrE35cQAvH17B+bfCYx2e48Y7GnkdSaRYU+EWEZE8ERERz4ZWlQG4cO0UUnRJbk/s3PU1faaOg9Aswtf35aVR93sdSaTYU+EWEZE80+TuvmSGQOsdu5k8epXXcYodn+8ICzteySuf7SBmax1m9B+BVjUX8Z4Kt4iI5Jk6De/k69MjCXOQ9N5gr+MUO7OGt6fTwn30+BX+GfEvLr6gpNeRRAQVbhERyUPh4WVIvrw686vDwuiV7N3rdaLiY8uKDzj3mW8AeLphax4e0dvjRCLyBxVuERHJUxc8OIjEPuFMbbGS/4xL8jpOsZCVeYDV3W+nUqpjfqVornpvhlb7FClAVLhFRCRP1azcg3pWDczxxrejvY5TLMx8IJFLV6RwKAIW3vwxZ50V6XUkETmKCreIiOSp0NBIbrugAY23w83rXuK3VX6vIxVpa394mUveXArA8CbXcP/QyzxOJCLHUuEWEZE81/v8AXz6gTF40W4mPfmW13GKrIyM3XSf8gITzoJpp5Xj9vHjCNHf7CIFjn4sRUQkz5Urn8j8hDgAqn77In4d5A6KRyZey5KSm7itXQzrBn2v1SRFCigVbhERyXNmRpW7sk9t6LZ1NV/M2OdxoqJn6TfDeGPZDwA02vY8991ex+NEInIiKtwiIhIUF3R8kkVVwonJgIUjHvM6TpGSlrqZrF5P8PPbRzhncSs+f/pmzLxOJSInosItIiJBERlZg98CS703WzmeQ4c8DlREOOeYedMlJCRlEJURwgNd3qRSJbVtkYJMhVtERILm3AF3kh4KrXfsYcKbi72OUyT8PP1+rpi8CYAXLv4n199c1+NEIvJXVLhFRCRo6p59F580LM2bCfDugglexyn0Dh9aS+jtrxPpg3fOqMVjHz7jdSQRyQEVbhERCZqwsGi2PH4Rd3WE7yrPZO1a53WkQss5P59d35rG2zPYVDqE0kPnEBfndSoRyQkVbhERCaobmw8kIrM0xP/KkFHfeR2n0Ppx1j/pND0JgFcTB3FNz9M9TiQiOaXCLSIiQVU+riXto6tww1KoN+MuMjK8TlT4pKSupevXH3N9N3j2rAY8+t4TXkcSkVOgwi0iIkFlZjyc0I53psLAlcuY9O4aryMVKs75uXPsdeyMTGJCjcqc/sy3xMZ6nUpEToUKt4iIBN25LYfwVe1IIvyw5T/3ex2nUPlx4m0sXboUgMsyXuWqK8p6nEhETpUKt4iIBF14eFn2da8PQLffZ7N2tc/jRIXDwT2/ENNvNIvezqTj122YOLyb15FE5G9Q4RYRkXyReM9zbI4xzjiYycQhI7yOU+A552Nur8s4c1cWm2LCuOO+sZQu7XUqEfk7VLhFRCRflK/Qmq8vLA9A/W9HkJ7ucaAC7vuxN9H5i534gbfaPcEVnap4HUlE/iYVbhERyRdmRv2Bvcky6Lw1iQ9fXep1pAJrf/Iiyv9zHKEOXmlwFoPfecjrSCKSCyrcIiKSb5q0eJy3m8VwZwcY9v1HOK2D8z/8/izm9exA3T1ZrIoLo8ELcylVyutUIpIbKtwiIpJvwsJKEz7sMt5KgA01R/P1t7oo97E+/6Ifrb5Lxmcw+vJnaNc+3utIIpJLKtwiIpKvel74b2LT46H0Dh4c/b7XcQqUvQeWcN28aTS5De5uegmD37rP60gikgdUuEVEJF9FRZ3BnbVq8OwseHHGPWze5Pc6UoHg92dy4zs3sj9yOxusFlc9O4OoKK9TiUheUOEWEZF81++Kh+m53Gi+M5UPB77odZwC4euXulHti+Xgh+4l/kObFtFeRxKRPKLCLSIi+a5ylc7MbVEBgIbznyU11eNAHtuxbhr1/jWT1z9z9J3djtGPt/U6kojkIRVuERHJd2ZGo8F3cSQMrti2k7cHz/A6kmd8vlSWX9WbKil+vqsUxS0vTSYy0utUIpKXVLhFRMQT55w/kGnnxQIQP/GeYrsQzqyHEmm77BAp4TC314ec31TXABQpalS4RUTEEyEhEZz2eC8yQ+CaLesZ++y3XkfKd+t+fIWLXl0EwLAmXXjkmS4eJxKRYFDhFhERzzRt8wzTG0YT6mDzxw/i83mdKP+kH9nB7h4PUiYdplUvy23jPiY01OtUIhIMKtwiIuKZsLBoYgZ3pfk/4ImOS3n3o11eR8oXzjmGv9ON8PQ0kqOMpP6fcfrpYV7HEpEgUeEWERFPtezwEjsq14LwIzw87blisdz7glXPMHz7LzS7GfpdNpg7+jf1OpKIBJEKt4iIeCo8vCzDLrsYgBqRrzFp1GqPEwXX/j0/0/W9l/GHHyZkU1deH/kvzLxOJSLBpMItIiKeu/qC5xn2TRl+GpPK/qd7FNlzuX2+I/zQvg0jp2wjbkc1xvUaRfnyatsiRZ0nhdvMyprZF2a2NnAbd5wxjc3sezNbYWbLzKy7F1lFRCT4IiIqkNCnNX7ghvVLef+FH7yOFBSzHmxJ+0X7ab8ObvcN5cr2Zb2OJCL5wKsj3AOBuc65OsDcwP1jHQZucM6dBbQHRphZmXzMKCIi+ahNj1FMOyeaCD9EvHYDmZleJ8pbK+cOpcWrPwEw5NxODH39Jo8TiUh+8apwdwHGBLbHAF2PHeCcW+OcWxvY3gYkAxXyLaGIiOSr8PAyxA2/lowQ6L5pLW//a6bXkfJMyt4VZPUeTukMmFAznrsmfkKYLkoiUmx4VbjjnXPbAQK3FU822MyaAhHA7/mQTUREPHLxFS8zOSGWEKD2mJvYv6/wX7LE709nfqeWNNyRwe+xofiGzqd6dV1wW6Q4MRek6y+Z2Ryg0nGeegQY45wrc9TYfc65/zmPO/BcZWAe0Mc5d9yT+sysL9AXID4+/rzx48fnMv3fk5KSQnR0tCf7lvyjeS76NMfe2r9jLIk3vkuZdHi49QDaPdYhKPvJr3ne/Mnd3PDyr6SHwt3th9BrQMug71P+pJ/n4sGreW7VqtVi51zCX40L2i+0nHOXnug5M9tpZpWdc9sDhTr5BONigBnAoycq24F9jQRGAiQkJLjExMRcZf+75s2bh1f7lvyjeS76NMfecu4SRk2aRdaGbbx12mdcX+kpGtTP+yPC+THPXy17kr5R6znQFJLozGufDCY8PKi7lGPo57l4KOjz7NUpJdOAPoHtPsDUYweYWQTwCTDWOTcxH7OJiIiHzELo+Ob73N+uHLtrreDqJ18rlIvhbNv9FV0+fIX0yMM8UK8L/adMUdkWKaa8KtxPAW3NbC3QNnAfM0sws7cDY64FLgFuNLOlga/G3sQVEZH8VKlcKx5qch4A68sN4YORv3qc6NSkpyYxtdvVhLIddtdnRt/3qFxZ19sWKa48+Yy0c24P0OY4jy8Cbglsvw+8n8/RRESkgBjYeQwL5pzLKx9vZ0tEG/Zds4O4sgW/tPqyjvDlZedyx3d7abQ+hC8GTqDNxaW9jiUiHtJKkyIiUiCVKFGJ4TfcStkj0GZHMqOuvsPrSH/JOT8zb2zC5d/tIi0UJlw4jMH9zvE6loh4TIVbREQKrPMuGsInvWoDcMuCkcx4b4XHiU5uzrDL6PDBagAGXtSLZ8c97HEiESkIVLhFRKTAMjN6vTSTL2pHUibdUeqBRJK3Z3kd67h+GnsHFz0xhxBgeKPzGfrpB/qQpIgAKtwiIlLAlSpVl5hRD7CzFCTu3M1HV3QucFctWfLlcM7o+x9KZcLoM6rS65MFxMR4nUpECgoVbhERKfCaXvI4n91zLgB3/vIZzw76xONEf1q+4X0unvUyo5s4PqlRlvrvr6RWLa3bLiJ/UuEWEZECz8zoPfRL3mxdnt7d4KHDDzL5s91ex2LN1ik0f3Mgh6OSGdCoCbGj13DhBTq0LSL/nwq3iIgUCmFhsVw5cSqz6lSFcuu4dnJXfl2V7lmeNXOfY2vz7kT6t2K76zGp2yxaJ5bzLI+IFFwq3CIiUmjEl72IL295nIi0sjQM+Y6lnRqwa1tGvudY/sljlOn8IK02ZTD0szjebzuXru0q5HsOESkcVLhFRKRQaVz7Zsa1782kj4zev2/gu2YN8MzEyQAACG1JREFU2JucmW/7n/9Md2r2eIKKhx2zTitN9aeX0Ktz1Xzbv4gUPircIiJS6Fx5yQgWPNaOAyWg6+b1zG/agOStaUHdp9+Xycwbm9B84ARKZ8C40ysQ+d5aruhQI6j7FZHCT4VbREQKHTOj1/2fMXVIGw5GQNdNv7MqoTprl+wMyv4OpqzjqyaVuWLMUkIdPNXwHBJmb6Fly/ig7E9EihYVbhERKZTMjOsf+oLpT3dieylouWMX/lY1mf7i9Dzdz5ylr1B9aFu+KbuHgxFwf+LV3PnNL9SpXSJP9yMiRZcKt4iIFFpmxnX9p/HLe/1ZXiGUMw6m8dxPvbj8gZGkpORudZzdq7+i/81NaDv5nxwotZFhjarz4oAPeP7LicTEWB69AxEpDlS4RUSk0Gt/5YuUXDCJm6+qxvx6KXwefRvl77uEsfc8Q1a675Rea9+a7/i0dX1KntOaR8YvJS49k7IbrmNe76UMHt4LU9cWkVOkwi0iIkXCGWd04a0Pl3Nb5asISStDi/RvueGVh9hYMZoPWrbj18nzcD7/cf/b1A3L+eL+7syrU4HS9S+m41erKZUJP1SMoX/UhySPep8W58fl8zsSkaJCa8+KiEiRER5ehjf7fswD275i5D/7syF2OWccSOOM+V/A/C84FG6sKxPNlrgY7m13DllThpGZsZUdb6ymbeA1MkNgUq0yLL70bh58ZgidYnVsSkRyR4VbRESKnNpVWvH0uF/YvO0Hnnv8Map8u5iLk/ZR/aCjya5DlHCH2Fh+63/H/1YONkSX4qfTzqJCz0HccmtXrgr38A2ISJGiwi0iIkVW9SoXMuA/XwBwMHUb4z4dT9KiFew+eIjE1Fjiy1bhzMr12DenOe0b1eBynZ8tIkGgwi0iIsVCTKkq9Ox+P3TPvj9v3jwSExM9zSQixYNOTBMRERERCSIVbhERERGRIFLhFhEREREJIhVuEREREZEgUuEWEREREQkiFW4RERERkSBS4RYRERERCSIVbhERERGRIFLhFhEREREJIhVuEREREZEgMuec1xnylJntAjZ5tPvywG6P9i35R/Nc9GmOiwfNc/GgeS4evJrnGs65Cn81qMgVbi+Z2SLnXILXOSS4NM9Fn+a4eNA8Fw+a5+KhoM+zTikREREREQkiFW4RERERkSBS4c5bI70OIPlC81z0aY6LB81z8aB5Lh4K9DzrHG4RERERkSDSEW4RERERkSBS4c4DZtbezFab2TozG+h1Hsl7ZnaamX1lZqvMbIWZ3et1JgkeMws1syVm9qnXWSQ4zKyMmX1sZr8Ffq6beZ1J8paZ3Rf48/pXMxtnZpFeZ5LcM7N3zCzZzH496rGyZvaFma0N3MZ5mfF4VLhzycxCgdeAy4EzgZ5mdqa3qSQIsoB/OucaABcCd2mei7R7gVVeh5Cgegn43DlXH2iE5rtIMbOqwD1AgnPubCAU6OFtKskjo4H2xzw2EJjrnKsDzA3cL1BUuHOvKbDOObfeOZcBjAe6eJxJ8phzbrtz7ufA9iGy/3Ku6m0qCQYzqwZ0AN72OosEh5nFAJcAowCccxnOuf3eppIgCANKmlkYEAVs8ziP5AHn3Hxg7zEPdwHGBLbHAF3zNVQOqHDnXlVgy1H3k1ARK9LMrCbQBPjR2yQSJCOABwG/10EkaE4HdgHvBk4detvMSnkdSvKOc24r8BywGdgOHHDOzfY2lQRRvHNuO2QfIAMqepznf6hw554d5zFd+qWIMrNoYBLQ3zl30Os8krfMrCOQ7Jxb7HUWCaow4FzgDedcEyCVAvgraPn7AufwdgFqAVWAUmbW29tUUpypcOdeEnDaUferoV9bFUlmFk522f7AOTfZ6zwSFM2Bzma2kezTw1qb2fveRpIgSAKSnHN//JbqY7ILuBQdlwIbnHO7nHOZwGTgIo8zSfDsNLPKAIHbZI/z/A8V7tz7CahjZrXMLILsD2VM8ziT5DEzM7LP91zlnHvB6zwSHM65Qc65as65mmT/LH/pnNNRsSLGObcD2GJm9QIPtQFWehhJ8t5m4EIziwr8+d0GfTC2KJsG9Als9wGmepjluMK8DlDYOeeyzKwfMIvsT0G/45xb4XEsyXvNgeuB5Wa2NPDYw865mR5mEpG/727gg8CBkvXATR7nkTzknPvRzD4Gfib7KlNLKOArEUrOmNk4IBEob2ZJwGDgKWCCmd1M9j+2rvEu4fFppUkRERERkSDSKSUiIiIiIkGkwi0iIiIiEkQq3CIiIiIiQaTCLSIiIiISRCrcIiIiIiJBpMsCiogUYmZWDpgbuFsJ8JG9bDnAYeecFvsQEfGYLgsoIlJEmNkQIMU595zXWURE5E86pUREpIgys5TAbaKZfW1mE8xsjZk9ZWbXmdlCM1tuZrUD4yqY2SQz+ynw1fwvXr+ymc03s6Vm9quZtciP9yUiUtjolBIRkeKhEdAA2Ev2yopvO+eamtm9ZK+62B94CXjROfetmVUnewXdBid5zV7ALOfccDMLBaKC+g5ERAopFW4RkeLhJ+fcdgAz+x2YHXh8OdAqsH0pcKaZ/fHfxJhZaefcoRO9JvCOmYUDU5xzS4MTXUSkcNMpJSIixUP6Udv+o+77+fPgSwjQzDnXOPBV9SRlG+fcfOASYCvwnpndEITcIiKFngq3iIj8YTbQ7487ZtY4cNvUzMYeO9jMagDJzrm3gFHAufkVVESkMFHhFhGRP9wDJJjZMjNbCdweeLw6cOQ44xOBpWa2BLiK7HPARUTkGLosoIiInJSZPQu855xb5nUWEZHCSIVbRERERCSIdEqJiIiIiEgQqXCLiIiIiASRCreIiIiISBCpcIuIiIiIBJEKt4iIiIhIEKlwi4iIiIgEkQq3iIiIiEgQ/R9YijbIgPY7JwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 864x360 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# First derivative of Gaussian function\n",
    "\n",
    "# Initiation of numerical and analytical derivatives \n",
    "nder_for=np.zeros(nt)      # forward FD operator\n",
    "nder_back=np.zeros(nt)     # backward FD operator\n",
    "nder_cent=np.zeros(nt)     # central FD operator\n",
    "ader=np.zeros(nt)          # analytical derivative\n",
    "\n",
    "# Numerical FD derivative of the Gaussian function\n",
    "for it in range (1, nt-1):\n",
    "    \n",
    "    nder_for[it]=(f[it+1]-f[it])/dt        # forward operator\n",
    "    nder_back[it]=(f[it]-f[it-1])/dt       # backward operator\n",
    "    nder_cent[it]=(f[it+1]-f[it-1])/(2*dt) # central operator\n",
    "\n",
    "# Analytical derivative of the Gaussian function\n",
    "ader=(-(time-t0)/a)*(1/np.sqrt(2*np.pi*a))*np.exp(-(time-t0)**2/(2*a))     \n",
    "\n",
    "# Plot of the first derivative and analytical derivative\n",
    "plt.plot(time, nder_for,label=\"FD forward operator\", lw=2, color=\"y\")\n",
    "plt.plot(time, nder_back,label=\"FD backward operator\", lw=2, color=\"b\")\n",
    "plt.plot(time, nder_cent,label=\"FD central operator\", lw=2, color=\"g\")\n",
    "plt.plot(time, ader, label=\"Analytical derivative\", ls=\"--\",lw=2, color=\"red\")\n",
    "plt.title('First derivative')\n",
    "plt.xlabel('Time, s')\n",
    "plt.ylabel('Amplitude')\n",
    "plt.legend()\n",
    "plt.grid()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The approximation of the first derivative of the Gaussian by all three FD operators seem to be very accurate. How large are actually the errors?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "code_folding": []
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAtwAAAFNCAYAAAAze7gSAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJzs3XlYVGX7wPHvA6jgkkqpmahYucuqiBu4o2WSqaRmKlkuqb2V/Vyy3FtMK8vybU/NPe11yTQ1hVyyBBPNXSvcFxRFEVCW5/fHgZFlWFSGI3B/rutcM3PmnOfc5xxG73nmWZTWGiGEEEIIIYRt2JkdgBBCCCGEEEWZJNxCCCGEEELYkCTcQgghhBBC2JAk3EIIIYQQQtiQJNxCCCGEEELYkCTcQgghhBBC2JAk3EIIkQOlVA2lVKxSyr6AjjdJKbXgLvZfp5QakJ8xpSs7Vin1sC3KFkKIokwSbiGEAJRSkUqp+NSkMm15SGt9QmtdVmudfAdlBiulttki3uxorR/TWs+723KUUqFKqRcylV1Wa/3P3ZYthBDFjSTcQghxS9fUpDJtOZPTxspwT/w7ei/FYitKKYe8rMtDOQXya4UQQqQp0v84CyHE3VJKuSqldFpil1rz+7ZSajsQBzycWpP9j1LqmlLqX6VUX6VUfeBzoHlqbfmVbMqvpZT6NXXfjcADmd5vppT6TSl1RSm1RynVJt171mIJVUq9oJQqlbpPo3TbV0qtxa+slKqolFqjlIpSSl1Ofe6Sut3bgB/waWrsn6au10qpR1NjOpc+cVVKPaWU2pv63E4pNVYp9bdS6pJS6nullHMO1/gJpVREary/KaXc070XqZQak1r2daWUQzbr6qee+xWl1H6lVGC6MuYqpT5TSq1VSl0H2iqlHldKHUi97qeVUv+X4x+CEELcBUm4hRDi9vUDBgPlgChgFvCY1roc0AKI0FofBIYCO1JryytkU9YiYBdGoj0VsLS/VkpVA34C3gKcgf8DflBKVcomluNpK7XWN4D/AX3Sbfs08KvW+gLGv/9zgJpADSAe+DR13zeArcCI1NhHpA9Ya/07cB1ol271M6nnAvAfoBvQGngIuAzMtnbySilv4FtgCHA/8AWwWilVKt1mfYAuQAWtdVLmdYACfgQ2AJWBl4CFSqm6meJ7O/U6bQO+AYak3rNGwGZr8QkhRH6QhFsIIW5ZmVpDekUptTKH7eZqrfenJn9JQArQSCnlpLU+q7Xen5eDKaVqAD7AeK31Da31FozEMc2zwFqt9VqtdYrWeiMQDjxuLRatdWKmQywiY8JtSYq11pe01j9oreO01tcwktHWeYk71eK0spVS5VJjWpz63hDgDa31qdTEfxLQM5vmH4OAL7TWf2itk1Pbn98AmqXbZpbW+qTWOj6bdc2AssA0rfVNrfVmYE2mc1+ltd6eeh0TgESggVLqPq31Za31n7dx7kIIcVsk4RZCiFu6aa0rpC7dctjuZNoTrfV1oBdGbfZZpdRPSql6eTzeQ8Dl1DLSHE/3vCYQlO5LwBWgFVDVWixWbAaclFK+SqmagCewAkApVVop9YVS6rhS6iqwBahwG+2bFwHdU2uiuwN/aq3TYq8JrEgX80EgGahipZyawGuZzrE6xrXJ6RzTr3sIOKm1Tkm37jhQLYcyemB8STie2qSneU4nK4QQd0MSbiGEuH06wwut12utO2IkwoeAr6xtZ8VZoKJSqky6dTXSPT8JzE/3JaCC1rqM1npadrFkiisF+B6jpvcZYE1qbTbAa0BdwFdrfR/gn7pe5SV2rfUBjKT2MTI2J0mL+7FMcTtqrU9bKeok8HambUtrrRen28ZaLOnXnQGqZ+o0WgM4nc32aK3DtNZPYjRBWYlxnYQQwiYk4RZCiLuglKqilApMTZpvALEYtbkA5wEXpVRJa/um1giHA5OVUiWVUq2Aruk2WQB0VUp1UkrZK6UclVJt0jo35tEijBr4vmRMisthtNu+ktqhcWKm/c4DuY25vQijvbY/sCzd+s+Bt1Nr1dM6az6ZTRlfAUNTa+GVUqqMUqpLajOVvPoDo035aKVUidSOpV2BJdY2Tr3WfZVS5VOb4Vzl1j0TQoh8Jwm3EELcHTuM2uIzQDRGO+hhqe9tBvYD55RSF7PZ/xnAN3XficB3aW9orU8CTwLjMDpnngRGcRv/dmut05LRh4B16d76CHACLgK/Az9n2vVjjHbXl5VSs7IpfjHQBtistU5/fh8Dq4ENSqlrqeX7ZhNfOEY77k8xOlceA4LzeHppZdwEAjFq2y8C/wX6a60P5bBbPyAytTnNUIz28kIIYRNK69x+8RRCCCGEEELcKanhFkIIIYQQwoYk4RZCCCGEEMKGJOEWQgghhBDChiThFkIIIYQQwoYk4RZCCCGEEMKGrE2zW6g98MAD2tXV1ZRjX79+nTJlyuS+oSjU5D4XfXKPiwe5z8WD3Ofiwaz7vGvXrota60q5bVfkEm5XV1fCw8NNOXZoaCht2rQx5dii4Mh9LvrkHhcPcp+LB7nPxYNZ91kpdTwv20mTEiGEEEIIIWxIEm4hhBBCCCFsSBJuIYQQQgghbKjIteEWQgghRN4kJiZy6tQpEhISzA7FZsqXL8/BgwfNDkPYmK3vs6OjIy4uLpQoUeKO9peEWwghhCimTp06Rbly5XB1dUUpZXY4NnHt2jXKlStndhjCxmx5n7XWXLp0iVOnTlGrVq07KkOalAghhBDFVEJCAvfff3+RTbaFyA9KKe6///67+iVIEm4hhBCiGJNkW4jc3e3nRBJuIYQQQpjG3t4eT09PyxIZGUloaCjly5fHy8uLunXr4u/vz5o1a6zuf+PGDTp06ICnpydLly4t4OizatOmjWnzgaQXGRnJokWLzA7jtq1cuZIpU6YAsGXLFry9vXFwcGD58uUZtps3bx61a9emdu3azJs3L8cyb968ib+/P0lJSTaLOzfShlsIIYQQpnFyciIiIiLDusjISPz8/CxJdkREBN26dcPJyYn27dtn2Hb37t0kJiZmKSMnycnJ2Nvb33XsSUlJODiYl0rldPy0hPuZZ57Jl/Jykvl65vX6Wjve9OnTWb16NQA1atRg7ty5vP/++xm2iY6OZvLkyYSHh6OUonHjxrRt2zbbNtwlS5akffv2LF26lL59+97u6eULqeEWQhRvWsNff8G8efD66zBkCPXeeQeGD4fJk2HlSjh71uwohSjWPD09mTBhAp9++mmG9RcuXODZZ58lIiICT09P/v77bzZt2oSXlxdubm4MHDiQGzduAMZM1FOmTKFVq1YsXLiQxo0bA7Bnzx6UUpw4cQKARx55hLi4OH788Ud8fX3x8vKiQ4cOnD9/HoBJkyYxePBgAgIC6N+/P/Hx8fTu3Rt3d3d69epFfHy81XPIKa4xY8bQtGlTmjZtyrFjxwCIioqiR48e+Pj44OPjw/bt260eP+3Libe3N97e3vz2228AjB07lq1bt+Lp6cnMmTNJSEjgueeew83NDS8vL0JCQgCYO3cuQUFBdO3alYCAgCxxL1iwgKZNm+Lp6cmQIUNITk4GoGzZskyYMAFfX1927NiR4fouW7aMiIgImjVrhru7O0899RSXL18GjF8Axo0bR+vWrfn4448zHOvIkSOUKlWKBx54wHJt3N3dsbPLmK6uX7+ejh074uzsTMWKFenYsSO//PILx48fp3bt2ly8eJGUlBT8/PzYsGEDAN26dWPhwoXZ/YnZnNRwCyGKH61h1y749lurCfWDABs3ZtynWTPo0weeew5kxAMh8k18fDyenp4A1KpVixUrVljdztvbmxkzZmRYV7lyZb7++mvef/991qxZQ0JCAm3atGHTpk3UqVOH/v378/XXXzN27FjAGNpt27ZtALz33ntcvXqVrVu30qRJE7Zu3UqrVq2oXLkypUuXplWrVvz+++8opfj666+ZPn06H3zwAQC7du1i27ZtODk58eGHH1K6dGn27t3L3r178fb2zhJ7QkICwcHBGeL67LPPeOWVVwC477772LlzJ9999x2vvPIKa9as4eWXX+bVV1+lVatWnDhxgk6dOlmGvUt//Li4ODZu3IijoyNHjx6lT58+hIeHM23aNMt1ASyx//XXXxw6dIiAgACOHDkCwI4dO9i7dy/Ozs4Z4j548CBLly5l+/btlChRgmHDhrFw4UL69+/P9evXadSokaX5R+br6+7uzieffELr1q2ZMGECkydP5qOPPgLgypUr/Prrr1mu0/bt261ev8xOnz5N9erVLa9dXFw4e/YsNWvWZMyYMQwdOhRfX18aNGhg+RLRqFEjwsLCci3bViThFkIUH1rDhg0wZQqk1gIBULUqtGoFDRrAgw9yMDKS+i4ucPIk7N4N27fD778by8SJMHIk/N//gZOTeeciRD4LDbVN58k2bXSO71trUmKN1jmXA3D48GFq1apFnTp1ABgwYECGWtRevXpZnrdo0YLt27ezZcsWxo0bx88//4zWGj8/P8AYMrFXr16cPXuWmzdvZhgOLjAwEKfUz/+WLVv4z3/+AxhJpru7e57imj17tiXh7tOnj+Xx1VdfBeCXX37hwIEDljKuXr3KtWvXshw/MTGRESNGEBERgb29vSWJzmzbtm289NJLANSrV4+aNWtatk2rLc5s06ZN7Nq1Cx8fH8D4clS5cmXAaHvfo0ePDNunXd+YmBiuXLlC69atLecbFBSUZbvMzp49S6VKlay+l561v4W0To0vvPACy5Yt4/PPP8/wd2Vvb0/JkiVNGyZSEm4hRPFw5AiMGHGr5trZGfr3NxZPT0jXA/18aCj127S5te/16/DTTzBrlpF8T5gAc+bAZ59Bp04Fex5CFFO7d++mfv36OW6TW1JepkwZy3M/Pz+2bt3K8ePHefLJJ3nvvfdQSvHEE08A8NJLLzFy5EgCAwMJDQ1l0qRJVsuB3EewyC2u9PunPU9JSWHHjh2WxDq785g5cyZVqlRhz549pKSk4OjoeNsxZD6f9PsMGDCAd999N8t7jo6OWdppZ1dOXo/n5ORETExMrvu7uLgQGhpqeX3q1Cl8fX0BiIuL49SpUwDExsZmSK5v3LiR7fWxNWnDLYQo2pKSjBptNzcj2a5YEd57D44fh5kzwcsrQ7JtVZky8PTTsG0bhIQYZf37L3TubNR037xZMOcihA21aaNtsuSHvXv3MnXqVIYPH57jdvXq1SMyMtLSDnr+/Pm0bNnS6rb+/v4sWLCA2rVrY2dnh7OzM2vXrrVsHxMTQ7Vq1QByHAXD39/f0jZ437597N27N09xpdX+ApbRVZYuXUrz5s0BCAgIyNBmPbtfAWJiYqhatSp2dnbMnz/f0sa6XLlylhrxzHEeOXKEEydOULdu3WzPC6B9+/YsX76cCxcuAEZnxePHj+e4DxizPlasWJGtW7daPd/s1K9f33KNctKpUyc2bNjA5cuXuXz5Mhs2bLB0ph0zZgx9+/ZlypQpDBo0yLLPpUuXqFSp0h3PFHm3JOEWQhRdJ09C27ZGM5CbN43214cPw+jRULbsnZXZpg38+Se88w7Y28MHHxi13KkdgoQQ+WPr1q2WYQGHDx/OrFmzsoxQkpmjoyNz5swhKCgINzc37OzseP75561u6+rqChiJKECrVq2oUKECFStWBIzOiUFBQfj5+Vk68Vnz4osvEhsbi7u7O9OnT6dp06Z5imvo0KGW92/cuIGvry8ff/wxM2fOBGDWrFmEh4fj7u5OgwYN+Pzzz60ef9iwYcybN49mzZpx5MgRS+2xu7s7Dg4OeHh4MHPmTIYNG0ZycjJubm706tWLuXPnUqpUqRyvZ4MGDXjrrbcICAjA3d2djh07cjaPncjnzZvHqFGjcHd3JyIiggkTJuS6j7+/P7t377bUxoeFheHi4sKyZcsYMmQIDRs2BMDZ2Znx48dbOpROmDABZ2dnfv31V8LCwixJd8mSJZkzZw4AISEhPP7443mK3RZUXtpEFSZNmjTRZo1/GRoaSpv0P0OLIknucyHxyy/QqxdERxtttBcsgHbt8rRrnu/xjh3Qo4fR6bJ+fVi/HtJ15BH3NvksG53icmumUdjd61O7u7q6Eh4enmNSX5y8/PLLdO3alQ4dOtzWfrnd5+7du/Puu+/mWqufE2ufF6XULq11k9z2lRpuIUTR8+WXRnOP6Gh47DHYsyfPyfZtad7c6EjZsCEcPGjUfqe2HRRCCHH7xo0bR1xcXL6WefPmTbp163ZXyfbdkoRbCFF0aA1jx8KQIZCcDGPGwJo1kIde73esRg2jbXfjxvDPP0YTlnPnbHc8IUSREhkZKbXb6VSpUoXAwMB8LbNkyZL0798/X8u8XZJwCyGKhpQUeOklo0OkgwN88w1MmwZ2BfDPXIUKRodMb284dgy6djVGNhFCCCGQhFsIURSkpBi12rNnQ8mSsGIFDBxYsDFUrAjr1kGtWhAeDs88Y8QlhBCi2JOEWwhRuKWkwPPPw9dfg6Mj/PgjpI6jW+AqVzaS7ooVYfVqYyQTIYQQxZ4k3EKIwktreO01mDsXSpeGtWshdRpf09StC4sWGWN7T5gAmzaZG48QQgjTScIthCi83n0XPvoISpSAlSuNDov3gs6d4c03jS8EffrA6dNmRyTEPcve3h5PT0/LEhkZSWhoKOXLl7eMw+3v78+aNWus7j9p0iTef//9u44jODiY5cuX33U5eeXq6srFixcL7HjZiYiIYO3atWaHcds++ugjvvvuOwCWLVtG06ZNsbOzI/PQ0O+++y6PPvoodevWZf369TmWGRUVRefOnW0Sr0ztLoQonL76Ct54w6hJXrAAOnY0O6KMJk6E334zariff95oapLbjJZCFENOTk5ZZlGMjIzEz8/PkmRHRETQrVs3nJyccp385l6UlJSEg4N5KVdOx4+IiCA8PPy2JoW5k/PRWqO1xi5dR/bk5OQs08Nbk3m7pKQkvv32W/78808AGjVqxMKFCxk5cmSG/Q4cOMCSJUvYv38/Z86coUOHDhw5ciTbY1aqVImqVauyffv2bGcovVNSwy2EKHw2bIAXXzSe//e/xrTr9xp7e+OLgLOzMSHON9+YHZEQhZanpycTJkzIMNV5env27KFdu3bUrl2br776CoDY2Fjat2+Pn58fbm5urFq1yrL9d999h7u7Ox4eHvTr1y9LeePHjyc4OJidO3fSvXt3AFatWoWTkxM3b94kISGBhx9+GICvvvoKHx8fPDw86NGjh2UM6eDgYEaOHEnbtm0ZM2YMly5dIiAgAC8vL4YMGUJ2Ew8uXrwYNzc3GjVqxJgxYyzry5Yty2uvvYa3tzft27cnKioKgL///pvOnTvTuHFj/Pz8OHTokNXj79y5kxYtWuDl5UWLFi04fPgwN2/eZMKECSxduhRPT0+WLl1KdHQ03bp1w93dnWbNmlmmqp80aRKDBw8mICDA6hB7M2bMwMfHB3d3dyZOnAgYX5zq16/PsGHD8Pb25uTJk5QtW5YJEybg6+vLjh072LRpE15eXri5uTFw4EBu3LgBGL8ATJkyhVatWrFs2bIMx9q8eTPe3t6WpL9+/frUrl07S0yrVq2id+/elCpVilq1avHoo4+yc+dOwsLCcHd3JyEhgevXr9OwYUP27dsHQLdu3Vi4cKHVe3NX0r5xFJWlcePG2iwhISGmHVsUHLnPJjt4UOvy5bUGrceNs8kh8vUeL1xoxFqunNbHj+dfueKuyWdZ6wMHDpgdgrazs9MeHh7aw8NDd+vWTWtt3JsuXbpk2G737t26Xr16WfafOHGidnd313FxcToqKkq7uLjo06dP68TERB0TE6OvXr2qo6Ki9COPPKJTUlL0vn37dJ06dXRUVJTWWutLly5prbUeMGCAXrZsmR41apQePHiwTklJ0YmJidrV1VVrrfVrr72mmzRpordt26ZDQ0N17969tdZaX7x40RLLG2+8oWfNmmUpr0uXLjopKUlrrfVLL72kJ0+erLXWes2aNRqwxJDm9OnTunr16vrChQs6MTFRt23bVq9YsUJrrTWgFyxYoLXWevLkyXr48OFaa63btWunjxw5orXW+vfff9dt27a1evyYmBidmJiotdZ648aNunv37lprrefMmWMpS2utR4wYoSdNmqS11nrTpk3aw8PDcp29vb11XFxclnuwfv16PWjQIJ2SkqKTk5N1ly5d9K+//qr//fdfrZTSO3bssGwL6KVLl2qttY6Pj9cuLi768OHDWmut+/Xrp2fOnKm11rpmzZr6vffey3IsrbWeMGGC5TqnuXr1qm7durUOCwuzrBs+fLieP3++5fXAgQP1smXLtNbGvXrttdf0sGHD9DvvvGPZ5tSpU7pRo0ZWj2vt8wKE6zzkp9KkRAhReERHG2Ncx8RA9+4wdarZEeWuTx9YvtwYqvCFF4zabmlaIu5BtvqzzKYi18JakxLr5WRf0JNPPomTkxNOTk60bduWnTt30qVLF8aNG0doaCgODg6cPn2a8+fPs3nzZnr27GmZbMbZ2dlSztSpU/H19eXLL78EwMHBgUcffZSDBw+yc+dORo4cyZYtW0hOTsbPzw+Affv28eabb3LlyhViY2Pp1KmTpbygoCBL84UtW7bwv//9D4AuXbpQsWLFLOcRFhZGmzZtqJQ6WVffvn3ZsmUL3bp1w87Ojl69egHw7LPP0r17d2JjY/ntt98ICgqylJFWQ5z5+DExMQwYMICjR4+ilCIxMdHqtdy2bRs//PADAO3atePSpUvExMQAEBgYiJOTU5Z9NmzYwIYNG/Dy8gKMXxeOHj1KjRo1qFmzJs2aNbNsa29vT48ePQA4fPgwtWrVok6dOgAMGDCA2bNn88orrwBYzjezs2fPZpli3RprfzMq9Q99woQJ+Pj44OjoyKxZsyzvV65cmTNnzuRa9u2ShFsIUTgkJkLPnsbEMl5e8N13BTOpzd1SCj77DH791ZgcZ/lySPefoxAib3bv3p1tkqUyfVtQSrFw4UKioqLYsmULzs7OuLq6kpCQgNY6y/ZpfHx82LVrF9HR0ZZE3M/Pj3Xr1lGiRAk6dOhAcHAwycnJlo6awcHBrFy5Eg8PD+bOnUtoaKilvDJlyuQYZ2Y5famwds4pKSlUqFAh2y8s6Y8/fvx42rZty4oVK4iMjKRNmzZ5jiEt7sznk36f119/nSFDhmRYHxkZmWUfR0dHy5eA3M43u+M5OTmRkJCQ474ALi4unDx50vL61KlTPPTQQwBER0cTGxtLYmIiCQkJlmMlJCRY/VJxt0z930op1VkpdVgpdUwpNdbK+yOVUgeUUnuVUpuUUjXNiFMIcQ8YORJCQuDBB40xrrP5h/ieVKXKrTG5R46E2Fhz4xHCCqPtU/4v+WHv3r1MnTqV4cOHW31/1apVJCQkcOnSJUJDQ/Hx8SEmJobKlStTokQJQkJCOH78OADt27fn+++/59KlS4CReKXp3LkzY8eOpUuXLly7dg0Af39/PvroI5o3b06lSpW4dOkShw4domHDhgBcu3aNqlWrkpiYmGPbX39/f8v769at4/Lly1m28fX15ddff+XixYskJyezePFiWrduDUBKSoplFJVFixbRqlUr7rvvPmrVqmVp46y1Zs+ePVaPHxMTQ7Vq1QCYO3euZX25cuUs55o5ztDQUB544AHuu+++bM8LoFOnTnz77bfEpv7bdvr0aS5cuJDjPgD16tUjMjKSY8eOATB//nzL+eakfv36ln1yEhgYyJIlS7hx4wb//vsvR48epWnTpgAMHjyYqVOn0rdv3wxt5Y8cOUKjRo1yLft2mVbDrZSyB2YDHYFTQJhSarXW+kC6zXYDTbTWcUqpF4HpgPXfF4QQRdeiRfDpp8YskitXgouL2RHdvhdeMEZW2bUL3n7bGNJQCJGtrVu34uXlRVxcHJUrV2bWrFnZjlDStGlTunTpwokTJxg/fjwPPfQQffv2pWvXrrRu3Rpvb2/q1asHQMOGDXnjjTdo3bo19vb2eHl5ZUhAg4KCuHbtGoGBgaxduxZfX1/Onz+Pv78/AO7u7lSuXNlS65vWDKVmzZq4ubllSF7TmzhxIn369MHb25vWrVtTo0aNLNtUrVqVd999l7Zt26K15vHHH+fJJ58EjNre/fv307hxY8qXL8/SpUsBWLhwIS+++CJvvfUWiYmJ9O7dGw8Pjyxljx49mgEDBvDhhx/Srl07y/q2bdsybdo0PD09ef3115k0aRLPPfcc7u7ulC5dmnnz5uV2qwgICODgwYM0b94cMDp4LliwINcRSBwdHZkzZw5BQUEkJSXh4+PD0KFDcz3eY489lqGz64oVKxgxYgQXL16kS5cueHp6sn79eho2bMjTTz9NgwYNcHBwYPbs2djb2/Pdd9/h4ODAM888Q3JyMi1atGDz5s20a9eOkJAQunTpkmsMt0vdzs8X+XpgpZoDk7TWnVJfvw6gtbb6v5BSygv4VGud4zgtTZo00ZnHYCwooaGh2f5EI4oOuc8FbP9+aNoU4uKMphl5+Mf4btnsHv/xBzRrZowbvm8fpLZbFOaQzzIcPHgwT21hC7Nr165Rrlw5s8O4a2XLlrXUIAt46qmnmD59umV0kvy6z/7+/qxatcpqG3trnxel1C6tdZPcyjWzSUk14GS616dS12XneWCdTSMSQtxbrl2DHj2MZPvZZyFT+8BCx9fXGJM7MRFef93saIQQotCaNm0aZ8+ezdcyo6KiGDlypNVk+26ZWcMdBHTSWr+Q+rof0FRr/ZKVbZ8FRgCttdY3rLw/GBgMUKVKlcZLliyxaezZiY2NpWzZsqYcWxQcuc8FRGsaTJlC5dBQrru6suu//yXFBh1ZrLHlPS556RK+zz6LfUICf86ezdUGDWxyHJE7+SxD+fLlefTRR80Ow6byOrmKKNwK4j4fO3bMMmJLmrZt2+aphtvMUUpOAdXTvXYBsozDopTqALxBNsk2gNb6S+BLMJqUmPUTofw8WTzIfS4gn38OoaFQrhxlfv4Z/7p1C+zQNr/Hu3fD22/jvWSJMXqJDBNoCvksGz+RF4XmFjkpKk1KRM4K4j47Ojpahj68XWY2KQkDaiulaimlSgK9gdXpN0htt/0FEKi1zr27qxAc5ryRAAAgAElEQVSiaDhwAF591Xj+5ZdQgMl2gRg1Cu6/H7ZuhZ9+MjsaIYQQNmZawq21TsJoJrIeOAh8r7Xer5SaopQKTN1sBlAWWKaUilBKrc6mOCFEUZGQAL17G4/BwcbzoqZ8eXjzTeP5uHGQkmJuPEIIIWzK1IlvtNZrgbWZ1k1I97xDgQclhDDXmDHw119QuzZ88onZ0djOiy/CBx8Y57pypTFzphBCiCKpEEzTJoQoNn76CWbNMobNW7QIinKHtlKlYGzqfF9TpuTfDCFCFDL29vZ4enpalsjISEJDQylfvjxeXl7UrVsXf39/1qxZY9M4QkND+e233+5ovyeeeMIGEd2+uXPn2mRaclvSWtOuXTuuXr0KwMCBA6lcuXKWyWeio6Pp2LEjtWvXpmPHjlYnDkrv008/Zc6cOTaL+3ZJwi2EuDecPWs0IQFjYpgmuXb6Lvyefx4eegj27DFmzxSiGHJyciIiIsKyuLq6AsaU6rt37+bw4cPMmjWLESNGsGnTJpvFkVPCnZSUZLPj3q7k5ORs37uThDun8rKT+Xrk9fpYO9batWvx8PCwzGYZHBzMzz//nGW7adOm0b59e44ePUr79u2ZNm1ajscaOHAgs2bNylNcBUESbiGE+VJSjGT74kXo0AFee83siAqGo6PRhAaklluIHHh6ejJhwgQ+/fTTLO/Fxsby3HPP4ebmhru7Oz/88AMAGzZsoHnz5vj5+REUFGSZNMbV1ZWJEyfi7e2Nm5sbhw4dIjIyks8//5yZM2fi6enJ1q1bCQ4OZuTIkbRt25YxY8awc+dOWrRogZeXFy1atODw4cM5xpyQkGCJy8vLi5CQEMBIip988kk6d+5M3bp1mTx5smWfBQsW0LRpUzw9PRkyZIglQS1btiwTJkzA19eXHTt2MGXKFHx8fGjUqBGDBw9Ga83y5csJDw+nb9++eHp6Eh8fz6ZNm/Dy8sLNzY2BAwdy48YNyzWYMmUKrVq1skwLnyYqKooePXrg4+ODj48P27dvB2DSpEkMHjyYgIAA+vfvz9y5cwkKCqJr164EBASgtWbUqFE0atQINzc3y0yYoaGhtG3blmeeeQY3N7cs12nhwoWW2TTBmHjG2dk5y3arVq1iwIABAAwYMICVK1cC8J///IcpU6YAsH79evz9/UlJSaF06dK4urqyc+fOHO9TgdFaF6mlcePG2iwhISGmHVsUHLnPNvDpp1qD1g88oPWZM2ZHU7D3OC5O6wcfNM5/7dqCO66Qz7LW+sCBA2aHoO3s7LSHh4f28PDQ3bp101ob96ZLly4Zttu9e7euV69elv1Hjx6tX375Zcvr6OhoHRUVpf38/HRsbKy+evWqnjZtmp48ebLWWuuaNWvqWbNmaa21nj17tn7++ee11lpPnDhRz5gxw1LOgAEDdJcuXXRSUpLWWuuYmBidmJiotdZ648aNunv37tnGqrXW77//vg4ODtZaa33w4EFdvXp1HR8fr+fMmaMffPBBffHiRR0XF6cbNmyow8LC9IEDB/QTTzyhb968qbXW+sUXX9Tz5s3TWmsN6KVLl1rKvnTpkuX5s88+q1evXq211rp169Y6LCxMa611fHy8dnFx0YcPH9Zaa92vXz89c+ZMyzV47733rNwNrfv06aO3bt2qtdb6+PHjlms+ceJE7e3trePi4rTWWs+ZM0dXq1bNEsvy5ct1hw4ddFJSkj537pyuXr26PnPmjA4JCdGlS5fW//zzj9Xj1ahRQ1+9ejXDun///Vc3bNgww7ry5ctneF2hQgWttdbXr1/XDRo00GvWrNF16tTRx44ds2zz1ltv6ffff9/qce+Etc8LEK7zkJ+a2mlSCCE4dgxGjzaef/EFVK1qbjwFzcnJGAJxzBijE+Vjj5kdkSim1GTbjAevJ+b8y01ak5Jcy8nmF6BffvmF9BPeVaxYkTVr1nDgwAFatmxJSkoKSUlJNG/e3LJN99ROyo0bN+Z///tftscMCgqyTKYSExPDgAEDOHr0KEopEhMTc4x327ZtvPSSMZdfvXr1qFmzJkeOHAGgY8eO3H///ZZYtm3bhoODA7t27cLHxweA+Ph4KleuDBjt3Hv06GEpOyQkhOnTpxMXF0d0dDQNGzaka9euGY5/+PBhatWqRZ06dQCjVnj27Nm88sorAPTq1Svb63ngwAHL66tXr3Lt2jUAAgMDcUo3AVnHjh0ttdHbtm2jT58+2NvbU6VKFVq3bk1YWBj33XcfTZs2pVatWlaPFx0dfVfjZ5cuXZqvvvoKf39/Zs6cySOPPGJ5r3Llyhw6dOiOy85PknALIcyTnGw0JYmLg759i+9IHYMHw9SpsGmTMSnOHU6sIERRtnv3burXr59lvdYalWnyKK01HTt2ZPHixVYnRClVqhRgJLI5tT8uU6aM5fn48eNp27YtK1asIDIyMtdJk7L7ggBkiVcphdaaAQMG8O6772bZ3tHR0ZL4JyQkMGzYMMLDw6levTqTJk0iISHhto4PGc8tvZSUFHbs2JEhsc5un/SvczpedscCcHBwICUlBTu7nFs5V6lShbNnz1K1alXOnj1r+TIC8Ndff+Hs7Jyl/XpCQoLV8zCDtOEWQphn5kzYvt3oOFiUhwDMTYUKMGiQ8fyDD8yNRRRbeqK2yZIf9u7dy9SpUxk+fHiW9wICAjK07b58+TLNmjVj+/btHDt2DIC4uDhL7XJ2ypUrZ6nJtSYmJoZq1aoBRjvs3Pj7+7Nw4UIAjhw5wokTJ6ibOonXxo0biY6OJj4+npUrV9KyZUvat2/P8uXLuXDBmOcvOjqa48ePZyk3Lbl+4IEHiI2NZfny5VbPoV69ekRGRlquwfz582ndunWucWe+nnn59SHtfJcuXUpycjJRUVFs2bKFpk2b5rpf3bp1+eeff3LdLjAwkHnz5gEwb948S7vv48eP88EHH7Bt2zbWrVvHH3/8YdnnyJEjWUY7MYsk3EIIcxw4cGvyl6+/hooVzY3HbC+/DPb2sHQpnDxpdjRCmG7r1q2WYQGHDx/OrFmzaN++fZbt3nzzTS5fvkyjRo3w8PAgJCSESpUqMXfuXPr06UPz5s1p1qxZrk0LunbtyooVKyydJjMbPXo0r7/+Oi1btszTyB7Dhg0jOTkZNzc3evXqxdy5cy01661ataJfv354enrSo0cPmjRpQoMGDXjrrbcICAjA3d2djh07cvbs2SzlVqhQgUGDBuHm5ka3bt0sTVDAGOFj6NCheHp6orVmzpw5BAUF4ebmhp2dHUOHDs017lmzZhEeHo67uzsNGjTg888/z3UfgKeeegp3d3c8PDxo164d06dP58EHH8x1vy5duhAaGmp5nXbPDh8+jIuLC9988w0AY8eOZePGjdSuXZuNGzcyduxYtNY8//zzvP/++1StWpVvvvmGF154wfKlZPv27XTocG9M6aJy+8mhsGnSpIkODw835dihoaG5/sQkCj+5z/kgKQmaN4fwcHjhBfjqK7MjysC0e9ynDyxZAv/3fzBjRsEfv5iRzzIcPHjQajONosRakxIzzZ07l/DwcKsjrhRHZ8+epX///mzcuPGuysl8n3fv3s2HH37I/Pnz7zZEC2ufF6XULq11ruPYSg23EKLgTZtmJNs1akgTivTShkP88ktInQRCCCGKsqpVqzJo0CDLxDf55eLFi0ydOjVfy7wbknALIQpWRASkjTs7Zw6kTnYgMCb7adPGSLbz0EZUCFH4BAcHS+12Jk8//bRl4pv80rFjR8skSvcCSbiFEAXnxg0YMMBoUjJiBLRrZ3ZE957UYcT47DOZCEcIIYoISbiFEAVnyhTYuxcefdRoViKyCgw0Rm05dAhSZ6YTQghRuEnCLYQoGH/8YSTZShnNJXIYl7VYc3CAIUOM5//9r7mxCCGEyBeScAshbC8+3mhKkpJidAxs2dLsiO5tgwYZiffKlXDqlNnRCCGEuEuScAshbO/NN+HwYahf35hRUeSsalVj1s3k5HtuyEQhbGHFihUope56Gu7g4OAME8FY884772R43aJFizs61qRJk3j//fdz3a5s2bK3Xfbjjz/OlStXbnu/K1eu8N90v4ydOXOGnj173nY5Iv9Jwi2EsK2tW40ZJe3t4bvvwNHR7IgKh2HDjMcvv4TERHNjEcLGFi9eTKtWrViyZInNj5U54f7tt99sfsy80lqTkpLC2rVrqVChwm3vnznhfuihh3L9AiIKhiTcQgjbiY2F4GBjtI1x44xh70Te+PtDgwZw7pzRtESIIio2Npbt27fzzTffZEi40yYm6tmzJ/Xq1aNv376kTdY3ZcoUfHx8aNSoEYMHDybzJH6bNm3iqaeesrzeuHEj3bt3Z+zYscTHx+Pp6Unfvn2BjDXQ06dPx83NDQ8PD8aOHQvAV199hY+PDx4eHvTo0YO4uLgcz+fff/+lefPm+Pj4MH78+AzvzZgxAx8fH9zd3Zk4cSIAkZGR1K9fn2HDhuHt7c3JkydxdXXl4sWLjBkzJkMCPWnSJD744ANiY2Np37493t7euLm5sWrVKsCYjfHvv//G09OTUaNGERkZaZna3NfXl/3791vKatOmDbt27eL69esMHDgQHx8fvLy8LGWJ/CUJtxDCdkaPhn/+AU/PW9O4i7xR6lbnydSpjYUoilauXEnnzp2pU6cOzs7O/Pnnn5b3du/ezUcffcSBAwf4559/2L59OwAjRowgLCyMffv2ER8fz5o1azKU2a5dOw4ePEhUVBQAc+bM4bnnnmPatGk4OTkRERHBwoULM+yzbt06Vq5cyR9//MGePXsYPXo0AN27dycsLIw9e/ZQv359y1Tj2Xn55Zd58cUXCQsLyzC1+YYNGzh69Cg7d+4kIiKCXbt2sWXLFgAOHz5M//792b17NzVr1rTs07t3b5YuXWp5/f333xMUFISjoyMrVqzgzz//JCQkhNdeew2tNdOmTeORRx4hIiKCGZlmq+3duzfff/89YMzueObMGRo3bszbb79Nu3btCAsLIyQkhFGjRnH9+vUcz1HcPkm4hRC2sXGjMZZ0iRIwbx6ULGl2RIVP377GdduwAU6eNDsaUdQpZZslF4sXL6Z3796AkRQuXrzY8l7Tpk1xcXHBzs4OT09PIiMjAQgJCcHX1xc3Nzc2b96coebWOBVFv379WLBgAVeuXGHHjh089thjOcbxyy+/8Nxzz1G6dGkAnJ2dAdi3bx9+fn64ubmxcOHCLMfKbPv27fTp0weAfv36WdZv2LCBDRs24OXlhbe3N4cOHeLo0aMA1KxZk2bNmmUpy8vLiwsXLnDmzBn27NlDxYoVqVGjBlprxo0bh7u7Ox06dOD06dOcP38+x7iefvppli1bBtxK3NPimjZtGp6enrRp04aEhAROnDiRY1ni9jmYHYAQogiKiYGBA43nkyaBu7up4RRa998PTz0FS5caQylm+nlaiMLu0qVLbN68mX379qGUIjk5GaUU06dPB6BUqVKWbe3t7UlKSiIhIYFhw4YRHh5O9erVmTRpEgkJCVnKfu655+jatSsAQUFBODjknPJorVFWviAEBwezcuVKPDw8mDt3LqGhobmel7VytNa8/vrrDEn75SpVZGQkZXIYJrVnz54sX76cc+fOWb6YLFy4kKioKHbt2kWJEiVwdXW1eg3Sq1atGvfffz979+5l6dKlfPHFF5a4fvjhB+rWrZvreYk7JzXcQoj89+qrxnB2TZsazUrEnUv74vLtt8awikLYita2WXKwfPly+vfvz/Hjx4mMjOTkyZPUqlWLbdu2ZbtPWmL5wAMPEBsbm22nwIceeoiHHnqIGTNmEBwcbFlfokQJEq10RA4ICODbb7+1tNGOjo4G4Nq1a1StWpXExMQszVCsadmypaUtevrtO3XqxLfffktsbCwAp0+f5sKFC7mW17t3b5YsWcLy5cstI47ExMRQuXJlSpQoQUhICMePHwegXLlyXLt2Lceypk+fTkxMDG5ubpa4PvnkE0s7+N27d+cak7h9knALIfLXjz/CnDlQqpTRlCSXWiWRiw4doEYNiIyUmSdFkbN48eIMnRsBevTowaJFi7Ldp0KFCgwaNAg3Nze6deuGj49Pttv27duXatWq0aBBA8u6wYMH4+7ubuk0maZz584EBgbSpEkTPD09LUP+TZ06FV9fXzp27Ei9evVyPaePP/6Y2bNn4+PjQ0xMjGV9QEAAzzzzDM2bN8fNzY2ePXvmmBynadiwIdeuXaNatWpUrVrVcl7h4eE0adKEhQsXWuK6//77admyJY0aNWLUqFFZyurZsydLlizh6aeftqwbP348iYmJuLu706hRoywdPUX+UJl79hZ2TZo00eHh4aYcO61HtSja5D7n4OJFaNQIzp+HDz80aroLoXvuHk+aBJMnwzPPQB5q2ETe3HP32QQHDx6kfv36ZodhMyNGjKB+/foMHz7c7FCEjV27do1y5crZ9BjWPi9KqV1a61yH4JIabiFE/hk+3Ei2W7eGl182O5qiIzjY6Hz2ww9w+bLZ0QhRKDRu3Ji9e/fSq1cvs0MRQhJuIUQ+WboUvv8eypQxmpTYyT8v+cbVFdq3hxs3oAAmBhGiKEgbdi99x0shzCL/Iwoh7t7Zs7dmRvzwQ6hVy9x4iqK0Tl8LFpgahhBCiNsnCbcQ4u5oDYMGQXQ0dO5sPBf5r1s349eD336Dv/82OxpRhBS1vlxC2MLdfk4k4RZC3J05c+Cnn6BCBfj66zxNdCHuQJkyxpjcIB0nRb5xdHTk0qVLknQLkQOtNZcuXcLR0fGOy5DxuoQQdy4yEl55xXj+6adQrZqp4RR5/foZTUoWLDAmwZEvN+Iuubi4cOrUKcsU6EVRQkLCXSVKonCw9X12dHTExcXljveXhFsIcWdSUoxJWa5dgx49jCHrhG21awcPPghHj0JYmDGxkBB3oUSJEtQq4n0uQkND8fLyMjsMYWP3+n2WJiVCiDsze7YxEUulSvDZZ1LbWhAcHKBPH+P5/PnmxiKEECLPJOEWQty+w4dhzBjj+ZdfGkm3KBjPPms8LlkCVqanFkIIce+RhFsIcXtu3jSaj8THQ//+xugZouB4eUGDBsasnhs2mB2NEEKIPJCEWwhxe8aPhz//NMba/uQTs6MpfpS6VcstY3ILIUShIAm3ECLvNm+GGTPA3t4Ymu6++8yOqHhK66C6cqXRaVUIIcQ9TRJuIUTeREcbTUi0Nmq5mzc3O6Liq2ZNaNUKEhLgxx/NjkYIIUQuJOEWQuROaxg8GE6fNhLtN94wOyLRq5fxuHSpuXEIIYTIlSTcQojczZ0LP/wA5coZ7YYdZAh/0/XsCXZ28PPPcOWK2dEIIYTIgSTcQoicHTsGL71kPJ89Gx5+2Nx4hOHBB6F1a2PUmFWrzI5GCCFEDiThFkJk78YN6N0brl83HtNGxxD3hqefNh6//97cOIQQQuRIEm4hRPZGjYJdu8DVVWaTvBf16GGMGLNhg9GpVQghxD3J1IRbKdVZKXVYKXVMKTXWyvv+Sqk/lVJJSqmeZsQoRLH1ww/GONslShgd8ypUMDsikVmlStCuHSQlwYoVZkcjhBAiG6Yl3Eope2A28BjQAOijlGqQabMTQDCwqGCjE6KY++cfGDjQeD5jBjRtam48IntpzUpktBIhhLhnmVnD3RQ4prX+R2t9E1gCPJl+A611pNZ6L5BiRoBCFEs3bhhJ3NWrxrTt//mP2RGJnHTvbowas3kzREWZHY0QQggrzEy4qwEn070+lbpOCGGm9O22v/1W2m3f65ydoWNHSE6G//3P7GiEEEJYYeZgutb+F9d3VJBSg4HBAFWqVCE0NPQuwrpzsbGxph1bFJyifJ8r/forDT/5hBQHB3aPHs21PXvMDskUhe0eV3F3p/66dVz+4gv21K1rdjiFRmG7z+LOyH0uHu71+6y0vqMc9+4PrFRzYJLWulPq69cBtNbvWtl2LrBGa708t3KbNGmiw8PD8znavAkNDaVNmzamHFsUnCJ7n/fvB19fYwjAjz6Cl182OyLTFLp7fPkyVK5szAh67hw88IDZERUKhe4+izsi97l4MOs+K6V2aa2b5LadmU1KwoDaSqlaSqmSQG9gtYnxCFF8XbkCTz1lJNt9+ki77cKmYkVjtJLkZFgt/4wKIcS9xrSEW2udBIwA1gMHge+11vuVUlOUUoEASikfpdQpIAj4Qim136x4hSiyUlKMCW2OHgUPD/j6a2m3XRj16GE8SjtuIYS455jZhhut9VpgbaZ1E9I9DwNcCjouIYqVyZPhp5+MzncrVkDp0mZHJO7Ek0/C0KGwcaMxwsx995kdkRBCiFQy06QQxdmqVTBlCtjZwZIlUKuW2RGJO1WlCvj5wc2bxhcoIYQQ9wxJuIUorvbuNZqSALzzjjG0nCjcunc3HqVZiRBC3FMk4RaiODp7Fp54AmJjjU6So0ebHZHID089ZTyuXQtxcebGIoQQwkISbiGKm7g4CAyEkyehRQuZ3KYoqVEDfHyMe7xhg9nRCCGESCUJtxDFSUoK9OsH4eFGe+2VK8HR0eyoRH6SZiVCCHHPkYRbiOJk3DgjEStfHtasgUqVzI5I5Le0hHv1aqMDpRBCCNNJwi1EcfHJJ/Dee2BvD8uXQ4MGZkckbKFOHWjUCGJiICTE7GiEEEIgCbcQxcPixbdmj/zyS+jQwdx4hG2lTYLzww/mxiGEEAKQhFuIom/9eujf33j+3nswcKC58QjbS2tWsnKlMd27EEIIU0nCLURR9vvvRvKVlASvvQajRpkdkSgIbm7w8MMQFWX8DQghhDCVJNxCFFX79kGXLsYQcQMGwIwZMvxfcaGUMdU7GLOJCiGEMJUk3EIURfv3Q7t2EB1tTHDz9deSbBc3knALIcQ9QxJuIYqagweNZDsqCgICYNkycHAwOypR0Fq2BGdnOHIEDh0yOxohhCjWJOEWoig5dAjatoULF6BjR5nYpjhzcDCaFIHUcgshhMkk4RaiqNi3z0i2z5+H9u2NZNvJyeyohJmkWYkQQtwTJOEWoij44w/w94dz54zmJKtXQ+nSZkclzNapE5QqZYxUcv682dEIIUSxJQm3EIXdpk1Gjfbly0aN5k8/SbItDGXLGn8bWsOPP5odjRBCFFuScAtRmP3vf/D443D9OvTrZ0zZLm22RXrSrEQIIUwnCbcQhZHW8MEH0LMn3LwJI0bA3LkyGonIqmtX4/GXX4wvZkIIIQqcJNxCFDZJSTBsGPzf/xmJ9zvvwKxZYCcfZ2FF1arg6wsJCbBxo9nRCCFEsST/QwtRmMTEQGAgfP650RluyRJ4/XWZ1EbkTJqVCCGEqSThFqKwOHAAmjaFdevggQdg82bo1cvsqERhEBhoPK5ZA8nJ5sYihBDFkCTcQhQGy5YZyfaRI+Dubgzz1qKF2VGJwqJBA3jkEbh4EX77zexohBCi2Mk14VZKlVZKjVdKfZX6urZS6gnbhyaE4OZNo632008bHd6eeQZ27DCSJyHySilpViKEECbKSw33HOAG0Dz19SngLZtFJIQwHDli1GJ/8IEx+sjHH8OCBTLGtrgz6RNurc2NRQghipm8JNyPaK2nA4kAWut4QHpoCWErWsPXX4OXF+zaBa6u8Ouv8J//SOdIcedatID774djx+DgQbOjEUKIYiUvCfdNpZQToAGUUo9g1HgLIfLb6dPQrRsMGgRxcdC3L0RESHttcfccHOCJ1NaA0qxECCEKVF4S7onAz0B1pdRCYBMw2qZRCVHcpKTAl18andtWr4b77oP5840mJOXLmx2dKCrSRiuRad6FEKJA5TotndZ6o1LqT6AZRlOSl7XWF20emRDFxYEDMHw4hIYar7t2hc8+g2rVTA1LFEEBAVCypDHKzfnzUKWK2REJIUSxkG0Nt1LKO20BagJngTNAjdR1Qoi7ceUKvPqqMcxfaChUqmRMZLNqlSTbwjbKloX27Y1+Aj/9ZHY0QghRbORUw/1B6qMj0ATYg1HD7Q78AbSybWhCFFGJifDttzB+PERFGR0hhwyBt94yJrQRwpYCA43Jk1avhoEDzY5GCCGKhWxruLXWbbXWbYHjgLfWuonWujHgBRwrqACFKDJSUmDRIqOd9tChRrLt5wd//mlM1S7JtigIaR0nN26E+HhzYxFCiGIiL50m62mt/0p7obXeB3jaLiQhipjkZGOmSE9PY9SRY8egdm1YutQY7s9TPk6iALm4QOPGxig4mzebHY0QQhQLeUm4DyqlvlZKtVFKtU6dcVIGcRUiN/HxRufHunWNmSL/+guqV4dvvjE6Sj79tIyrLcyRNlrJ6tXmxiGEEMVErqOUAM8BLwIvp77eAnxms4iEKOxOnjQmrvnsM6PZCMDDD8Nrr8Hzz0OpUubGJ4q8mzfh8uVbS3y8sS5tKePQlceYSPz3P7Ku42eUKGVHyZLGACblykHFisZSvjzY25t9NkIIUfjlZVjABGBm6iKEsCYpyeiI9sUXxmNKirG+SRMYPRq6d5fMRdw1reHcOaNV0okTcOqUMVfSqVPGcu4cREfD9eu5leTJCVyofuUU7wbtIhyfbLcsX96YoNLFJevy8MPw6KNQpky+nqYQQhQ5uSbcSql/SZ1lMj2t9cM2iUiIwkJrCAuDxYuN9thnzxrrS5QwmosMHQr+/tJsRNy2GzeM2df/+gsOH4ajR28tsbG5729vf6uWumJFKF3aqL0uVYrUmmzFofBAqh/9L6Prrmbeoz6W2u+rV2/VjMfE3Fr++Sf747m4QJ06xlK3rjHSpYeHkagLIYTIW5OSJumeOwJBgLNtwhHiHpeczH0HDsCGDUaSnT4LqV0bBg+GAQOMMbWFyIXWcOYM7N2bcTl0yPjRxBpnZ+NPzdXVSHSrVbv1WLWqkeSWK5eH73nrA6HzfwkqtZqgNVOtbpKcbCTbUVEZa9JPnTJq2P/+21jS1mXug1mtmpF4e3gYfYObNvi1yNoAACAASURBVIWaNeU7qBCi+MlLk5JLmVZ9pJTaBkywTUhC3GMuXzYS7J9+gnXr8L6YbqLVqlWhVy/o3dvIJiSTEDm4ehXCw+GPP24t585l3U4po6bYzQ3q1zcS7LQl32qN27QxJsLZuxeOHzcy4Uzs7Y0E39nZiMeapCRj9yNHjBr4/fthzx6jdv70aWNZu/bW9pUrg6/vrcXHx2i2IoQQRVlempSkn1XSDqPGu5zNIhLCbBcvwtatsGWLMWzfnj232mQD8Q8+iFOPHtCzpzGOdjFrm611MsnJ10lOvk5KSpzlefrXWt9E6yS0Trb6qJQdSpVAKQfLo51d2qMj9vblLIuDw32W53Z2eflR7t6QlAT79mVMrg8eNGq106tQwaj9dXe/tTRsaDQDsalSpaBzZ1i+HH78EUaMuKNiHBzgkUeM5bHHbq1PTjZqv/fsMZY//4SdO+HCBeNwP/5obKcU1KsHLVsaHyc/P6P2Xr67CiGKkrz87/VBuudJwL/A07YJR4gCdvUqRETArl1GRhAebvyen56Dg9EWu0sX6NKFP86do03btubEm0+SkxNISoomKekyiYmX8/DceExKikHrm6bFbW9/HyVLVqZEicqULFkl9dF47ehYA0dHVxwda+LgULBVplobg9OkT6537co6r0yJEkZynb6G99FHTUwuAwONhHv16jtOuLNjb3+rXXdQkLFOayMJ37nz1nXavdv4InLwoDG4DxhNUdKSbz8/4wuIXV4GsRVC/H97dx4nR1Xvffzzq+6efSY7WQghIUHCHkgAWdSwuYEsIhdkMaKAPsJFfFwQeK56uYiAeBFRr3ARDKggCkLYBFkCyBYChCUQIIQt+z5LZnp6O88fp3p6MpmZTJjpqVm+79frvM6p6uqu30xNT/+66tQ50kd1JeH+unNus9tlzGxSkeIRKY7aWn/32aJFhbJwob8O3lZZGRx4oE+yP/lJ+PjHNz/duGpV78XdiVwuHSbB+eTYJ8hdaedyyW7s2YjFKgmCSmKxinbaFQRBaXj2OhbWcSDWsg5y5HLp8Kz35nUu10QmU082WyiZTF3YrqOpqY6mps4nu43Hh1JWNpHS0h0pL9+JioqpYdmVkpLu96+vrd28a8i8ee13DZk8uZBY77+/T7bLyrq9+57z+c/7THbuXP9DFblvh5n/gjFlCpxyil/X3Oy/8z75pC//+pfvhnLbbb6Av/Gz9Rnw6dP9zZ/9iXNZ0ul1pFKrSKfXkMnUksnUks3WtmrXkcnUkcs1kcs1k8s141wzuVyyZRlyne7Hv8dKwvdgCUHQut16XQWxWBWxWGVYqlq9l7dcH4tVEQSV/eoqk0hf0pV3zt+AfdtZN727OzezzwLXADHgBufc5W0eLwVuDve1DjjJOfded/crA4xzPln48EN/J9cHH/hOpfn2kiWFEUTaSiR8R9np033Zd19/Tb8Xxsp2zpHLNYbJZV2bxHlDmzPLmz+WTq8nl9vq2G8dMksQjw8nkRhGPD6ceHxYJ21f+zKEICjDIjgl65wjk9kYJiyrSaVWt9Sp1Aqamz8gmXyPZPJ9MpmNNDQsoKFhwRavE4+PoKJiKpWVu1JZuQdVVftQVTWNeLym3f2mUr6b87x5hbJo0ZZdQ4YNKyTW+XrkyGL8JnrQiBE+k33ySXjwQT+6Ti8rLS18Kfne93zvrTfeKCTgTz7p39r33usL+C8tBxzgk+9DDvHfj2vaP3y9IpdrJpn8kGTyPZqb3yeZ9CWVWkkqtYpUaiXp9Bq2liz3B4UuX1Vhab8dj299G9+uwEyXL2Tg6zDhNrOpwO7AEDP7YquHavCjlXSL+dNcvwGOBJYCz5vZHOfc6602+zqwwTk3xcxOBq4ATuruvqWPyuX8AMJ1dT6BrqsrlPxQCWvW+E6gq1cX2mvW+KyoM2Vl/tr21Kmbl9133+qpMp8YJ8nlmshmG4GlNDS8TDbb2LLOn5FqJJtt3OKsbOuzs23Xde8DOCAeH7pZQtx+O59YF9pBUBFJ0twe5xw5lyPrsmRzWTK5DFkX1rksWZcl53KtSgk5N45cbAy5IEeuxK/H5ShxOWK5LKn0OpLNy2lKLiPZ/CGNTe/RlHyfpuT7ZLLrcDxFzj2FA3LOj3taWjKWsrIp1NftxUsvjeK/767l7bfLeWexkU4buACcAUZ8sjFlsrHbrsZuuxm7727sMN4IzIgHceJBnBXZOGvXJlqW40GcRNBmOeaXgygTjmOO8VntnDmRJNxtBYF/W+6+ux9ZE/z359YJ+Btv+NsrHn+88Jy99/bJ9yGH+ER87NiejSuXy5BMvktj4yIaG98M60Ukk0tIpVbSzui5W4jHh4ddoUYRjw8lHh9CPD6EWGxIq3YNsVj+KlEpQVAoZqXhFaKOuPAqUQrnUmHdvNly/qy5/1+1iWy2gVxuU0u7cE9G2/UN4X0aSXK5ZPgFomcEQeU2JOjtJfRVBEFZeDKgtKUdBCVK5qXPMNf2NE3+AbNjgeOAY4DW8//WA7c5557u1o7NDgR+4pz7TLh8IYBz7mettnkw3OYZ89ekVwKjXEdBAzNmzHDz58/vTmjbLJ1MsmLhqzz/wnNM2XkK2WwOcDjAnAOXBcDwn+wOh7kcOOf7bjqHcznM3ObrcjnM/GsYhEmF868DkH8dwPKJm3Pkfz1GFv+Szj/u/HN8EIXXaFl2jiCXhWwWy5dMFjJZLJfDspnCulyWIJuDTMY/lslBJkuQzUImTaw5jTWnCZpTWHOKoDntSzK/nPKPJ1MEyRSxhiZim5r97+sjyFQkSI6upmlMlS+jy2kcXU7jmHI2jSll03YJnGVwuTTOpcnm0kCaXC7flSFNzqXJ5dLgMuRcilw2STaXxLnmlv20F17bVfnlfCLnnE+r26uxBBaUh6WKIFaJBflS0VKCoAKLVWBWDuH2WMK/VpiIOlybxHTz0jqZ7Six7dLjnbTbe35XXmuwM2yzBDwRJCiJlVAaL6U0VkpZvKyl3WHdyWPtPT+/rvq9Fexy8DFkhw5h7bsLKSuratmmr3wpa2vtWnjqKd/95F//8t172g6juNNOhTPghxziR1lp++PMnTuXmTNnbrbOOUcy+R4NDS+zadPLNDS8TGPjIpqaFuNcuoOIAkpLx7fcQ5CvS0rGUlIympKSMSQSowiCftYPpg1/8qGpVQJeTzbbEF6l23xdR+222+ZyjUWM11/Ny7k48UQlZmXhl5cysJIwQS/DghJfLA7EwWKAv4mbfDc4ixOEyxBvecx3lUtA2H3OPzcWJvoBZgHOWbhsm62HIPyjDML3mmG0Xb/l8/zHtYWfR9byOZNf73/2Vo8712Y5/P1g5HLhYy2/NGv1WoQ5RX4HQA6cy+Fy5nMGXEsbcuGHW/i8nI8uv72FuU7h+DifKeVfOy8Xrg+3yR/MwmdrmPc4aPVyrFy7hlPPOHvb/1C6ycxecM7N2Np2HZ7hds7dDdxtZgc6557p0ei87YEPWy0vBQ7oaBvnXMbMaoERwFr6kOceuZ9Djj6BCVEHMgA0JKCu1JfaskK7rhTWVMCaSlgdljUVYV0JyUQaWB+WVlaFpc9Kh6WOPh5oURlGLIgRD+LELLZFO2YxAgs+UjECmpMByaaApkZfGuoDNjUE4Vnr8My1+S+fJaVpKipqqaxaT0XlRioqagli/kuBIyCIVRGL1RCEI6dg8ZYPDuccWZclnU2TyWVaSjrXZrnV4+lcGocjnUuTznWU0BXXGyNg6rpa/u0743liYmF9Saykw8S90/ZHeE5JrMQf8/DY549/e+tKq+N8+vNxPne0X59sCpg3zyffTz4JTz/te5ItWQKzZ/ufZeRI3w1lxgw/FOF++wGkqa9/gfr6l1qSa3/1qq7d31Np6Q7hvQC7UFExlfLyXSgt24lEyTiwoOWLbesvus0uS1MmRy69frPH819G29adPZa/2tPRY8V87Z6NK0E2V0PWVYaPZcmG74esy7R5Xj75coX3GbSc/SgkZi2nj9rI/49t2uKRyDgoy0BVCqqboTq1Zbsq5bcpzYR1trDcup1/LJGFeA5iLqxz27YcOP8v0PB1f7wucNvU7SGChLurOutS8gPn3JXAKWb25baPO+fO6+a+2zt10va90pVtMLOzgbMBRo8ezdy5c7sZ2rZZ+v4HTKoKv1WGITsrBOoAZ9aqTfvt/Padbevyr93+fgrrOooj37bNts8/ngmMrPk6E0DWjIz5tn+s7XK4rbXaPjCaYkYyHtAUN5IxIxkLSMaNprBOxgLfjhnJuFGXiFFfEpBtOxSBa+9PoNW6ZiC5lW06fK2ubNPOdlvZJggcRkBgYBYQmO9qELOAIDBiYYkHEMQgHhjxmBGL+XXxmJGIG/EYYe2Ix/1rBOGZkLa1YYVlLEw0C7WZ+aSVwCexYQLbksgSbLau0/Wt1nW0vr3ndLa+u10qUilj9eoyVq4sY8WKMpYtK+fDDyv44IMKli8vJ5fb8pgFgWOHHRqZPLmBKVN8mTy5geHD0zQ0NFBVVY4flOl1YGFYL8V/OWqdkG0P7NGqTGBbP66yLrtZyeQypF2aVC7Vkoh3uuzSLetTuZRf33qda/W8cJvWz3909zVMfWITJ75dyguT44XHsylS2RT1qfqPcFR6V8vf1CEB9gmjNGfkskYuZ2RzxtqccR9wX7PBv4AnDbMcsYdzxGIZgiBDLJYlMId/P8eAgBwBOZe/YrWanFtJjkf92d4B0Cd70HDW0iWsZbltu+V/+zZum4ORTY5xDY7t6x1jG3KMbMoxosn5OukY0Wp5WNIR/2gXcyOTo5BPtF9bJ48V8qDWeUhrW1wp7sI2bbdLlpb1ev63LTrrUvIF59w9Zjarvcedc7O7teMB1KUkr73Lk4Nd6yPV6spQi1zOl6zvydJhuzuPpdO+5Keu7qi9tceam6GxEVatqiWRGEJTkx/2rbGRlnauiJ+/lZV+npLq6p6pe+G+0G7L5WD9ej8CyMqVfoCYfHv5cnjvPV+WL+/4NYIAJk3avOv+3nvDHntAeXn7z+novZxKraW+/jlqa5+hru4Z6uqe2+Lm1Xh8GDU1BzJkyMEMGXIw1dX7EYsVe1DtbvrXv3z/i8mT/ew15i9Dp7IpkpkkzdlmmjPNNGeb/XI77eZM82bbbvV5bR5LZVMt3Y9ad0Fq3V2pvfVRXRVore0VmFjQZrmdKzT5Kzft1S1fSjvYprPnd/u53YgLApKNMRo3xWhsiLGpPkZDQ4yGuhgNdQH1dTHq62LUbYxRXxujdqNvb9wQo6kxgFwMXKxN7btUbHMybFBRbpSXQxAkGTq0jPJyWkpFRaFdVuZv5Skp8ffRt65LSiARdwxJrmLEhsUMW7+YoWsXU7P6HSo2LKV8/TJK1y0nlm5u+2fRKVdSQq6yGlfl/yG7quqWf85WXY1VVfrgSkuhvAwLa0pKsbJSrLzMP1YW1omEH8I2Hvdjcna1bl3M2i/9RFQ5WE90KbknrLuVWHfieWDncIjBZcDJwClttpkDzAKeAb4EPNpZsi19T+v3anvv2/44Z8zcuS+1+6Z2zifp7SXi+eWGBl/q6zevO2rn602bCqWnRiVMJDpPyCsr/f/x0lL/odNRXVKy+bFt7x2ayfgvLMnklnVTk78ntrYWNm70Jd+ure3al5hYDHbYwU+YMnGi77+bT6533rnnhuIrKRnJiBFHMWLEUYC/kW7TppeprX2qpaRSy1i//n7Wr/fTK5rFqaratyUBr6k5mNLSMT0TUE858EA/Ysk77/ghWHbdFTPzXT7iff+bWSabpLbuJTbUPcvG2nnU1c/fbOjI/J9kPD6C6uoZVFVNp6JiOg8+FMPsEyxYAC8tcCx6AzY1trom6ArdjrYbFWPyTgGTdwqYslOMnacETJkcMGmSMWLEwBkn3Dn/v2rDBl/Wry/U69ost603bmz//d8VsZgf6WfIED/iTHV1+6WzxyoqCol06/9Lc+c+27VEzDn/bf6VVwrltdf8l9BNWxkVauhQP4D89tvDuHG+D9OIEb60bo8YAcOHYyUl9MOPP+mGzrqU3EP7Z/ABcM4d050dh32yzwUexF+7u9E5t9DMLgHmO+fmAL8HbjGzxfjOuSd3Z58ixWRWSEB7ejjj/AAuW0vMt6VOpwsfqn3Z0KEwZkz7ZccdfYI9frw/YdPbgiBOdfV0qqunM378eTjnaG7+IEy+n6au7ikaGl6hvn4e9fXzWLr0agDKynbaLAGvrNwt2tEUYjE/sdPNN/vRSnbdNbpYtsLf1LiEurp51NU9R339c9TXv7TZzc0BUF1SRlXVvtTUHEB19f7U1BxAWdnEzW4EXTZ5LjNnDm/12n5ElNde81PT5+s334TVS3155oktY0okfI41blwh3xo3riW3YtiwQhk+3H+hLcaJQ+f8F9n6+s5LbW3HSfOGDVsf9KkzNTWFn7Mrdb5dXR3BydQNG+DZZ3155hk/+dm6de1vO3x4YRD5KVP81aAJEwoHvLKyd2OXfqezj6irir1z59z9wP1t1v2oVTsJnFjsOET6uiAonMXpqeHOUqnOE/JNmwpdaTqq8+222n5wBkHh6mhZqyuh+TJkiC9Dh25eDxnik5n+wszCUSp2ZPRof8Euk6mjru45amufoq7uKerqniWZXEIyuYRVq24B/EQ9+W4oNTUHU1Ozf+93QznmmELCfcEFvbvvTqTT66ire576+ueoq3uOurp5ZDJbJkUVFVOprj6AmhqfXFdW7kUQbNsfj1nhKsnRRxfWZ7N+Mp7Fi/1FgMWLC+XDD33e9v77vnRF/v1QVrZlnUgU3j9t63xC3fY9mC/J5JajtXwUpaVdT5hb10OHRvPFt8vWrYNHH4WHH4YnnthyVmHwP8Tee/v5Gfbay9e77OJ/SJFu6KxLyeP5tpmVAFPxZ7zfdFHO7SwiPaKkpHCFU4onHq9h+PAjGT78SCDfDeXVlgS8tvYpmps/ZP36B1i//gEg3w1ln5YEfMiQgykt7eGBpdv69Kf9H8Uzz/ix7Ud1f0bObZVOb2TTppepr3+J+vrnqa+f1+6soonEdpudua6u3o9EYmjR4orF/MnMCRPgsMO2fLypyd9HsGzZ5nXbM8f50thY6CLW0xKJjrtc5LuLDRnSeeLc0b0N/U42C08/zU7XX+9nVXrxxc37vJSW+mFrDjzQzyi8//7+clk/6rcs/cdWv4ua2VHA74B38HcmTDKzbzjnHih2cCIiA43vhrIP1dX7AOcCkEx+uFkC3tDwcphwPg/8EoCysknU1BwUnr3dk8rKPSgp6cGpLKurfTb5j3/AfffBV7/ac6/dhnM5kskPwqH4FrSUZPK9LbYNgjKqqqZTU3NASyktndCnxggvL/c9DCZP7tr26bQ/G52/h6F1nQ7vAW3vJnPY8r6K9sqg1twMjzwCf/+7v1qzenVhyN6SEj8o+xFH+L/1ffbZ6sRnIj2lKxd/fgEc6pxbDGBmk4H7ACXcIiI9oKxsB8rKTmb0aH+bSiZTH3afeCpMxJ8lmXyXZPJdVq/+U8vzEonRVFX55LuiYlfKyydTVjaZsrIdtjIjYQeOOcYn3HPm9EjCncnUkky+t9nMjPmZGtub8CQIyqis3JOqqmkt/a8rK/fY5q4hfV0iUTgTLT3AOd8Pe/Zs+Mtf/N2beTvtxIfTp7PDmWf6ZLuij48YJANWVxLu1flkO7QEWF2keEREBr14vJrhw49g+PAjAHAuS0PDq+FNmAtoaHiVxsaFpNOr2LBhFRs2PLzZ880SYV/yyZSWbk9JyZhw1sMxlJaOJZEYGU6TXU0sVlm4YfMLX4BvfQsefNCfbm0zvIufsS5NNttAOr2OdHoN6fQaUqk1pNNrSaWWkUy+31Ky2doOf8ZEYjuqqvYKk+t9qKqaRnn5xwiCvtwJWPqUVavghhvgD3/wHerz9toLvvhFOP542HNP3nn8cXbQkL0Ssa78Z1toZvcDt+P7cJ8IPG9mXwRwzt1ZxPhERAY9sxjV1dOorp7Wsq7QLeM1Nm16lcbGN0km36Gp6R1SqRU0NS1ut/9zO69OLFZFLFaFWZw9P5ag6q1GFv12EhsPqgSyZLON5HKNZLONsA2TvQRBOWVlO1Jevks4Q+PUlpkaEwndhCYf0XPPwbXXwu23F/rgjBkDp50GX/mKv9FRpI/pSsJdhp9z+lPh8hpgOPAFfAKuhFtEpJeZBZSXT6S8fCIjRx692WPZbCPJ5Ls0NS0hlVpBKrWypW5uXkEms45Mpp5stj5MpH0bYO2BUPUW1Mxdycp929tvgiCoIJEYQSIxipKSUSQSo0gkRlJSMrZllJbS0h1JJEb0qb7W0o/lcnDXXXDFFTBvnl8XBHDssfDNb/p+2X16iBQZ7Lb61+mcO6M3AhERkZ4Ri1VQWbk7lZW7b3Vb57Jksw1kMvVADkpegdlfYMz80Qzd73EsiBMElcRiFQRB+YDrTy19XDbrz2T/9KewcKFfN2wYnHmm7/40cWKk4Yl0VVdGKZkE/DswsfX23Z34RkREomcWIx4fQjweztZ04A6w/fYEy5ZRsagBpk+PNkAZnHI5fwPkT34Cb73l140f78eI/9rXdPOj9Dtduf5yF37Gx3vYls57IiLS/5j50Ur+53/8aCVKuKW3Pf64Hzd7/ny/PHEiXHghzJqlcQ+l3+rKXMJJ59yvnHOPOecez5eiRyYiItE4JryAOWdOtHHI4PLOO75P9syZPtkeOxauv96f4T77bCXb0q915Qz3NWb2Y+AhoDm/0jn3YtGiEhGR6Bx6qJ+ScMEC+OADP8WiSLE0N8PPf+77aSeTUFkJP/gBfPe7vi0yAHQl4d4TOB04jEKXEhcui4jIQFNa6qd6v/NOuOceOOecqCOSgWruXD/KyJtv+uVTT/XJ99ixkYYl0tO60qXkeGAn59ynnHOHhkXJtojIQKZuJVJMmzbBuef6qylvvgkf+5ifkv2Pf1SyLQNSVxLul4GhxQ5ERET6kM9/3o9z/NhjUFcXdTQykDz9NEybBr/5jR87+yc/gVdegcN0Lk8Grq4k3KOBRWb2oJnNCcvdxQ5MREQiNGoUHHSQn8nvoYeijkYGgkwG/uM/4BOf8FOx77GHn8Tmxz/WDZEy4HUl4f4xvlvJZcB/A/OAKcUMSkRE+gB1K5Gesny5nw3y0kv98gUX+JFI9tkn2rhEeslWE+5wCMBa4CjgD8DhwO+KG5aIiEQun3Dfd58/OynyUfzzn74LyeOPw5gx8PDDcPnlOqstg0qHCbeZfczMfmRmbwC/Bj4ELLxp8tpei1BERKKxyy7+Zrb1632/W5Ft4ZxPrD/zGVizBg4/3A81eeihUUcm0us6O8O9CH82+wvOuUPCJDvbO2GJiEif8IUv+FrdSmRbNDXBaaf5GSKd8zdGPvggjB4ddWQikegs4T4BWAk8Zmb/a2aHA9Y7YYmISJ+gftyyrZYvh099Cv78Zz9xzV13+RsjY7GoIxOJTIcJt3Pu7865k4CpwFzgO8BoM/sfM/t0L8UnIiJROuggGD4c3n67MDmJSEdeew323x+efx4mToRnnvHTtYsMcl25aXKTc+5PzrmjgfHAAuCHRY9MRESiF4/DUUf5ts5yS2cefxwOOQSWLfP1vHmw555RRyXSJ3RlWMAWzrn1zrnrNNOkiMggom4lsjV33OFvjqythRNO8COTjBoVdVQifcY2JdwiIjIIfeYzUFLiRypZsybqaKSvue46OPFEaG6Gc86Bv/wFysqijkqkT1HCLSIinauu9kO55XJw//1RRyN9yTXXwDe/6Uci+elP4dprdXOkSDuUcIuIyNZpeEBp68or4fzzffvXv4aLLgLTYGYi7VHCLSIiW5dPuB98EJLJaGOR6P3Xf/np2c3g+ut9VxIR6ZASbhER2boJE/z03Js2wdy5UUcjUbrsMvjRjyAI4Kab4Kyzoo5IpM9Twi0iIl2j0UrkV7+Ciy/2Z7ZvvhlmzYo6IpF+QQm3iIh0TeuE27loY5Hed+ON8O1v+/b118Opp0Ybj0g/ooRbRES6Zt99Ydw4P7HJSy9FHY30pttvL3QdufpqOPPMaOMR6WeUcIuISNeYqVvJYPTAA/5sdi4Hl1xSGJlERLpMCbeIiHSdhgccXF54wU9qk8nA978P/+//RR2RSL+khFtERLrusMOgosJ3KVm6NOpopJjeew+OOsqPTHP66XDFFRpnW+QjUsItIiJdV1bmp3oHuOeeaGOR4tmwAT73OVi1Cg4/HG64Qcm2SDco4RYRkW2jftwDW3MzHHccLFoEe+4Jd9wBJSVRRyXSrynhFhGRbXPUUf5s56OPQn191NFIT3IOzjgDnngCtt8e7r8fhgyJOiqRfk8Jt4iIbJtRo+CggyCV8lO9y8Bx+eVw661QVQX33Qfjx0cdkciAoIRbRES23fHH+/qOO6KNQ3rOPfcUZpH8859h772jjkhkwFDCLSIi2+6EE3x9zz3Q1BRtLNJ9CxfCKaf4LiWXXloY/lFEeoQSbhER2XYTJ8KMGX7IOHUr6d/Wr4djj4WGBjjpJLjwwqgjEhlwlHCLiMhH86Uv+fpvf4s2DvnoMhmfZL/zDuyzD9x4o4b/EymCSBJuMxtuZv80s7fDelgH2/3DzDaa2b29HaOIiGxFPuGeM8cPJSf9z49+BA8/DNttB3fd5Sc1EpEeF9UZ7h8CjzjndgYeCZfb83Pg9F6LSkREum7yZH9WtL4eHnoo6mhkW917L/zsZxAE8Je/wIQJUUckMmBFlXAfC8wO27OB49rbyDn3CKBBXkVE+ip1K+mf3n3XT9cOcNllMHNmpOGIDHRRJdyjnXMrAMJ6u4jiEBGR7sgn3Hff7cfllr4vmfTHbeNGPxrJ978fdUQiA168WC9sZg8DY9p56OIi7Ots4GyA0aNHM3fu3J7eRZc0NDREtm/pPTrOA5+O8baZsdNOVC1ZwitXX836Aw6IOpwuG6zHeeerr2b7yi7pCAAAEwBJREFUF1+kaexYXjjrLDJPPBF1SEU1WI/zYNPXj3PREm7n3BEdPWZmq8xsrHNuhZmNBVZ3c1/XA9cDzJgxw82M6NLY3LlziWrf0nt0nAc+HeNtNGsW/PjH7PXmm3DBBVFH02WD8jj/8Y/+JtfSUsrvvZdD9t036oiKblAe50Gorx/nqLqUzAFmhe1ZwN0RxSEiIt114om+vusuSKejjUU69uab8I1v+Pa118IgSLZF+oqoEu7LgSPN7G3gyHAZM5thZjfkNzKzJ4G/Aoeb2VIz+0wk0YqISMd23RV22w02bIDHHos6GmlPczOcfDI0NsKpp8KZZ0YdkcigEknC7Zxb55w73Dm3c1ivD9fPd86d2Wq7TzjnRjnnyp1z451zms5MRKQvyt88efvt0cYh7bvwQliwAHbaCX77W01uI9LLNNOkiIh037/9m6/vuEOT4PQ1DzwAV18N8TjceivU1EQdkcigo4RbRES6b/fdYa+9/FBzD+piZJ+xcqW/qRXgv/4L9t8/2nhEBikl3CIi0jO+/GVf33prtHGIl8v5ZHvNGjjsMPjBD6KOSGTQUsItIiI94+STfX333dDQEG0s4ruRPPQQjBgBt9zip3AXkUjo3SciIj1j4kQ46CBoavJjPUt0XnzR3ygJcNNNMG5ctPGIDHJKuEVEpOeoW0n0mprgtNP8mOjnnOOnbxeRSCnhFhGRnnPiiRCLwT/+AevWRR3N4HTxxfDGGzB1Kvz851FHIyIo4RYRkZ40ejQcfjhkMn6IQOldjz3m+27HYnDzzVBeHnVEIoISbhER6WnqVhKN2lr46ld9+z/+A/bbL9JwRKRACbeIiPSs44+H0lJ4/HFYtizqaAaP88+HDz6AGTPgoouijkZEWlHCLSIiPWvIEDj6aHAO/vSnqKMZHO66C/7wBygr80MAJhJRRyQirSjhFhGRnpef3XD2bJ94S/GsXg1nn+3bV1zhb5YUkT5FCbeIiPS8z34WRo2C11+H+fOjjmbgcs4n2/nZJM89N+qIRKQdSrhFRKTnJRJ+LGjwZ7mlOGbP9jN71tT4CW40m6RIn6R3poiIFEd+xIw//xmamyMNZUB6/3047zzf/vWvYcKEaOMRkQ4p4RYRkeLYay+YNg02bIB77ok6moHFOfj616G+3o8Kk7+aICJ9khJuEREpnvxZbnUr6VnXXQePPAIjR8LvfgdmUUckIp1Qwi0iIsVzyikQj8MDD8CqVVFHMzAsWQLf+55v//a3sN120cYjIlulhFtERIpn1Cg46ijIZjUmd0/I5eBrX4NNm+Ckk+DEE6OOSES6QAm3iIgUV75byU03aUzu7vrNb/wMnttt52+UFJF+QQm3iIgU1+c/7890v/YaPPts1NH0X2+/DRdc4NvXXef7b4tIv6CEW0REiqukBM44w7evuy7aWPqrbNb/Dpua/Igkxx0XdUQisg2UcIuISPGddZav//IXP0ygbJtrroGnnoKxY31bRPoVJdwiIlJ8U6bAEUdAMgm33BJ1NP3LokVw8cW+ff31MHx4tPGIyDZTwi0iIr3jG9/w9XXX6ebJrspk/E2nyaSvjz466ohE5CNQwi0iIr3j2GNh9Gh4/XXfPUK27he/gOeeg/Hj4eqro45GRD4iJdwiItI7Egk/hjTo5smuWLgQfvQj377hBhg6NNp4ROQjU8ItIiK956yz/DTkf/0rrF4ddTR9VzoNs2ZBKuV/Z5/5TNQRiUg3KOEWEZHeM2mSn3myuVlnuTtzxRXwwgswYQJcdVXU0YhINynhFhGR3nX++b7+7W/9GVzZ3MsvwyWX+PaNN0JNTbTxiEi3KeEWEZHeddhhsMcesHIl3H571NH0LamU70qSTsO3vgWHHx51RCLSA5Rwi4hI7zKDb3/bt6+5RkMEtnbppf4M96RJvluJiAwISrhFRKT3nXoqjBgB8+fD009HHU3f8NxzcNll/gvJTTdBVVXUEYlID1HCLSIiva+8HL75Td/+5S+jjaUv2LQJTjsNsln4v/8XPvWpqCMSkR6khFtERKLxrW9BPA533gnvvBN1NNH63vdg8WLYc0/46U+jjkZEepgSbhERica4cXD66ZDLwZVXRh1NdO6/H373OygpgT/+EUpLo45IRHqYEm4REYnOBRf4Pst/+AMsWxZ1NL1vzZrC7JuXXgp77RVtPCJSFEq4RUQkOrvsAiee6IfD+8Uvoo6mdzkHZ58Nq1bBJz/p+26LyICkhFtERKJ14YW+vu46WLs22lh60+zZcNddUF3t27FY1BGJSJEo4RYRkWhNm+ane29s9ONyDwZLlsB55/n2r38NEydGGo6IFJcSbhERid5FF/n62mthw4ZoYym2VApOPhnq6+GEE/yNoyIyoCnhFhGR6B10kJ/GvLYWfv7zqKMprosvhuefhx13hP/9X3/TqIgMaJEk3GY23Mz+aWZvh/WwdraZZmbPmNlCM3vFzE6KIlYREekl+fGnr7kGVq6MNpZieeABuOoq31/71lth2BYffyIyAEV1hvuHwCPOuZ2BR8LlthqBrzjndgc+C/zSzIb2YowiItKbDjgAjj3W9+W+9NKoo+l5y5fDV77i25deCgceGG08ItJrokq4jwVmh+3ZwHFtN3DOveWceztsLwdWA6N6LUIREel9l17qu1hcfz28+27U0fScbNZP3b52LRx5JPzgB1FHJCK9KKqEe7RzbgVAWG/X2cZmtj9QAgzyuX9FRAa4PfbwiWk6DT9s7+JnP3XJJfDYYzB6NNxyCwS6hUpkMDHnXHFe2OxhYEw7D10MzHbODW217QbnXLsd2cxsLDAXmOWce7aDbc4GzgYYPXr09Ntuu62b0X80DQ0NVFVVRbJv6T06zgOfjnG0SletYv9Zs4g1N/PS1VdTO21aUfbTW8d5xDPPsOdFF+GCgFeuvJIN06cXfZ9SoPfz4BDVcT700ENfcM7N2Np2RUu4O92p2ZvATOfcinxC7ZzbpZ3tavDJ9s+cc3/tymvPmDHDzZ8/v0fj7aq5c+cyc+bMSPYtvUfHeeDTMe4D/vM/4Sc/gb33hhdeKMqkML1ynN9+G/bbz4++cvnlfip76VV6Pw8OUR1nM+tSwh3VNa05wKywPQu4u+0GZlYC/B24uavJtoiIDBDf/z5MmAAvv+yHzuuPGhrg+ON9sn3CCeq3LTKIRZVwXw4caWZvA0eGy5jZDDO7Idzm34BPAl81swVhKc51RRER6VsqKvzweeDHrV69Otp4tpVz8PWvw8KFMHUq3HSTxtsWGcQiSbidc+ucc4c753YO6/Xh+vnOuTPD9h+dcwnn3LRWZUEU8YqISAS+9CU/osf69fDv/x51NNvmkkvg9tuhuhr+/ndfi8igpdukRUSkb8oPD1hZ6ZPXu+6KOqKu+dOffP/zIPCT20ydGnVEIhIxJdwiItJ3TZwIP/uZb3/rW7BhQ6ThbNWTT8LXvubbv/wlHHVUtPGISJ+ghFtERPq2c86Bgw6CFSvg7LN9/+i+aNEif5NkKgXnntv/usGISNEo4RYRkb4tCGD2bN8P+m9/891M+pr33oMjjoB16/xZ7auvjjoiEelDlHCLiEjfN2UKXHedb59/Prz6arTxtLZypb+5c9ky+MQnfH/zeDzqqESkD1HCLSIi/cOXv+yH2ksm4Ytf9GeTo7Z2LXz607B4Mey7L9xzjx/SUESkFSXcIiLSf/zqVzBtmk9wv/hFaG6OLpYVK+BTn/Jn26dOhX/8A4YMiS4eEemzlHCLiEj/UVHhzyKPGwdPPAFnnRXNTZTvvw+f/CS8/jrsvjs8+iiMGtX7cYhIv6CEW0RE+pfx4wtdN265Bb7znd5NuufPh49/vNCNZO5cGDu29/YvIv2OEm4REel/9t3Xj1hSUgLXXAPf/W7vJN133unPbK9cCTNnwiOPwMiRxd+viPRrSrhFRKR/+tzn4I47IJHww/Cddlrx+nSnUvC978EJJ0BTE5xxBjz4IAwdWpz9iciAooRbRET6r6OPhjlzoKoK/vxnOPxw+PDDnt3HokVw8MHwi19ALAZXXgm//70/uy4i0gVKuEVEpH/77Gf9lOrjxsFTT/lRTG6/vftdTBob4eKLYa+9fL/tHXf0+/n+98GsZ2IXkUFBCbeIiPR/06bBggW+m8n69XDSSX587Bde2PbXSibhN7+BXXaByy6DdNqP//3SS3DggT0fu4gMeEq4RURkYBg1Cu69F373Oxg+HB5+GGbM8N1Mbr4ZNm7s+LnptB9m8Lzz/Jnyc8+FpUt9Iv/UU3DDDTBsWO/9LCIyoGjuWRERGTiCAL7xDX9z4+WX++ngH33UFzPYdVfYbTcYOZKdly/3ifTixbBwITQ0FF5n+nS46CI47jj/miIi3aCEW0REBp6RI+Gqq3wf7Ntug7/+1fe/fv11X4Dt2z5n1119f/DTT4d99un1kEVk4FLCLSIiA9ewYfB//o8vySS88gq8+y6sXctbb73Fx/bZByZN8sn2dttFHa2IDFBKuEVEZHAoK4P99/cFWD53Lh+bOTPamERkUFDHNBERERGRIlLCLSIiIiJSREq4RURERESKSAm3iIiIiEgRKeEWERERESkiJdwiIiIiIkWkhFtEREREpIiUcIuIiIiIFJESbhERERGRIlLCLSIiIiJSROacizqGHmVma4D3I9r9SGBtRPuW3qPjPPDpGA8OOs6Dg47z4BDVcd7ROTdqaxsNuIQ7SmY23zk3I+o4pLh0nAc+HePBQcd5cNBxHhz6+nFWlxIRERERkSJSwi0iIiIiUkRKuHvW9VEHIL1Cx3ng0zEeHHScBwcd58GhTx9n9eEWERERESkineEWERERESkiJdw9wMw+a2ZvmtliM/th1PFIzzOzHczsMTN7w8wWmtm3o45JisfMYmb2kpndG3UsUhxmNtTM/mZmi8L39YFRxyQ9y8y+E/6/fs3MbjWzsqhjku4zsxvNbLWZvdZq3XAz+6eZvR3Ww6KMsT1KuLvJzGLAb4DPAbsBXzaz3aKNSoogA3zXObcr8HHgHB3nAe3bwBtRByFFdQ3wD+fcVGBvdLwHFDPbHjgPmOGc2wOIASdHG5X0kD8An22z7ofAI865nYFHwuU+RQl39+0PLHbOLXHOpYDbgGMjjkl6mHNuhXPuxbBdj/9w3j7aqKQYzGw8cBRwQ9SxSHGYWQ3wSeD3AM65lHNuY7RRSRHEgXIziwMVwPKI45Ee4Jx7AljfZvWxwOywPRs4rleD6gIl3N23PfBhq+WlKBEb0MxsIrAP8Fy0kUiR/BL4AZCLOhApmp2ANcBNYdehG8ysMuqgpOc455YBVwEfACuAWufcQ9FGJUU02jm3AvwJMmC7iOPZghLu7rN21mnolwHKzKqAO4DznXN1UccjPcvMjgZWO+deiDoWKao4sC/wP865fYBN9MFL0PLRhX14jwUmAeOASjM7LdqoZDBTwt19S4EdWi2PR5etBiQzS+CT7T855+6MOh4pioOBY8zsPXz3sMPM7I/RhiRFsBRY6pzLX6X6Gz4Bl4HjCOBd59wa51wauBM4KOKYpHhWmdlYgLBeHXE8W1DC3X3PAzub2SQzK8HflDEn4pikh5mZ4ft7vuGc+++o45HicM5d6Jwb75ybiH8vP+qc01mxAcY5txL40Mx2CVcdDrweYUjS8z4APm5mFeH/78PRjbED2RxgVtieBdwdYSztikcdQH/nnMuY2bnAg/i7oG90zi2MOCzpeQcDpwOvmtmCcN1Fzrn7I4xJRD66fwf+FJ4oWQKcEXE80oOcc8+Z2d+AF/GjTL1EH5+JULrGzG4FZgIjzWwp8GPgcuB2M/s6/svWidFF2D7NNCkiIiIiUkTqUiIiIiIiUkRKuEVEREREikgJt4iIiIhIESnhFhEREREpIiXcIiIiIiJFpGEBRUT6MTMbATwSLo4BsvhpywEanXOa7ENEJGIaFlBEZIAws58ADc65q6KORURECtSlRERkgDKzhrCeaWaPm9ntZvaWmV1uZqea2Twze9XMJofbjTKzO8zs+bAcvJXXH2tmT5jZAjN7zcw+0Rs/l4hIf6MuJSIig8PewK7AevzMijc45/Y3s2/jZ108H7gGuNo59y8zm4CfQXfXTl7zFOBB59xPzSwGVBT1JxAR6aeUcIuIDA7PO+dWAJjZO8BD4fpXgUPD9hHAbmaWf06NmVU75+o7ek3gRjNLAHc55xYUJ3QRkf5NXUpERAaH5lbtXKvlHIWTLwFwoHNuWli27yTZxjn3BPBJYBlwi5l9pQhxi4j0e0q4RUQk7yHg3PyCmU0L6/3N7Oa2G5vZjsBq59z/Ar8H9u2tQEVE+hMl3CIiknceMMPMXjGz14FvhusnAE3tbD8TWGBmLwEn4PuAi4hIGxoWUEREOmVmPwducc69EnUsIiL9kRJuEREREZEiUpcSEREREZEiUsItIiIiIlJESrhFRERERIpICbeIiIiISBEp4RYRERERKSIl3CIiIiIiRaSEW0RERESkiP4/sd5nHZOO0OEAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 864x360 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Erros of the FD approximation\n",
    "\n",
    "# Plot of the errors of the first derivative FD approximations\n",
    "plt.plot(time, 10*(nder_for-ader),label=\"FD forward operator error (10x)\", lw=2, color=\"y\")\n",
    "plt.plot(time, 10*(nder_back-ader),label=\"FD backward operator error (10x)\", lw=2, color=\"b\")\n",
    "plt.plot(time, 10*(nder_cent-ader),label=\"FD central operator error (10x)\", lw=2, color=\"g\")\n",
    "plt.plot(time, ader, label=\"Analytical derivative\",lw=2, color=\"red\")\n",
    "plt.title('First derivative errors')\n",
    "plt.xlabel('Time, s')\n",
    "plt.ylabel('Amplitude')\n",
    "plt.legend()\n",
    "plt.grid()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The errors for all three operators are actually quite small, so we have to amplify them by a factor 10. Larger errors occur in areas with significant slope. Notice also, that the central operator has the smallest error. You wonder why? We will see in a later lecture."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## FD approximation of the 2nd derivative\n",
    "\n",
    "We derived FD approximations for the 1st derivative. However, to approximate the acoustic wave equation, we also need to know how to approximate a 2nd derivative. \n",
    "To achieve this, we calculate the **Taylor expansion** of $f(t+dt)$ up to the second order term:\n",
    "\n",
    "\\begin{equation}\n",
    "f(t+dt) \\approx f(t) + \\frac{\\partial f(t)}{\\partial t} dt + \\frac{1}{2}\\frac{\\partial^2 f(t)}{\\partial t^2} dt^2 + \\mathcal{O}(dt^3)\\nonumber\n",
    "\\end{equation} \n",
    "\n",
    "We do the same for $f(t-dt)$:\n",
    "\n",
    "\\begin{equation}\n",
    "f(t-dt) \\approx f(t) - \\frac{\\partial f(t)}{\\partial t} dt + \\frac{1}{2}\\frac{\\partial^2 f(t)}{\\partial t^2} dt^2 + \\mathcal{O}(dt^3)\\nonumber\n",
    "\\end{equation} \n",
    "\n",
    "Adding the approximations $f(t+dt)$ and $f(t-dt)$ leads to \n",
    "\n",
    "\\begin{equation}\n",
    "f(t-dt) + f(t+dt) \\approx 2 f(t) + \\frac{\\partial^2 f(t)}{\\partial t^2} dt^2\\nonumber\n",
    "\\end{equation} \n",
    "\n",
    "Finally, we can rearrange to $\\frac{\\partial^2 f(t)}{\\partial t^2}$ and get the following **FD approximation for the 2nd derivative**\n",
    "\n",
    "\\begin{equation}\n",
    "\\frac{\\partial^2 f(t)}{\\partial t^2} \\approx \\frac{f(t-dt) + f(t+dt) - 2 f(t)}{dt^2}\\nonumber\n",
    "\\end{equation} \n",
    "\n",
    "##### Excercise\n",
    "\n",
    "* Calculate the 2nd derivative $\\frac{\\partial^2 f(t)}{\\partial t^2}$ of the Gaussian \n",
    "\\begin{equation} \n",
    "f(t)=\\dfrac{1}{\\sqrt{2 \\pi a}}e^{-\\dfrac{(t-t_0)^2}{2a}},\\nonumber\n",
    "\\end{equation}\n",
    "analytically. In case you are too lazy to compute the 2nd derivative by hand, you can also type **Second derivative of (1/sqrt(2*pi*a))*exp(-(t-t0)^2/(2*a)) with respect to t** in [Wolfram Alpha](https://www.wolframalpha.com/)\n",
    "* Compute and compare the numerical with the analytical 2nd derivative of the Gaussian and the errors between both solutions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "code_folding": []
   },
   "outputs": [],
   "source": [
    "# 2nd derivative of Gaussian function\n",
    "\n",
    "# Initiation of numerical and analytical derivatives \n",
    "nder2=np.zeros(nt)      # 2nd derivative FD operator\n",
    "ader2=np.zeros(nt)      # analytical 2nd derivative\n",
    "\n",
    "# Numerical FD derivative of the Gaussian function\n",
    "for it in range (1, nt-1):\n",
    "    \n",
    "    # ADD YOUR 2ND DERIVATIVE OF THE GAUSSIAN HERE!\n",
    "    # nder2[it]=        # 2nd derivative FD operator\n",
    "    \n",
    "# ADD ANALYTICAL 2ND DERIVATIVE OF THE GAUSSIAN HERE!\n",
    "# ader2=\n",
    "\n",
    "# Plot of the numerical and analytical second derivative of the Gaussian\n",
    "#plt.plot(time, nder2,label=\"FD operator\", lw=2, color=\"y\")\n",
    "#plt.plot(time, ader2, label=\"Analytical 2nd derivative\", ls=\"--\",lw=2, color=\"red\")\n",
    "plt.title('2nd derivative Gaussian')\n",
    "plt.xlabel('Time, s')\n",
    "plt.ylabel('Amplitude')\n",
    "plt.legend()\n",
    "plt.grid()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "code_folding": []
   },
   "outputs": [],
   "source": [
    "# Error of the 2nd derivative FD approximation\n",
    "\n",
    "# Plot of the second derivative FD approximation error\n",
    "#plt.plot(time, nder2-ader2,label=\"FD operator error\", lw=2, color=\"violet\")\n",
    "#plt.plot(time, ader2,label=\"Analytical 2nd derivative\", lw=2, color=\"red\")\n",
    "plt.title('2nd derivative errors')\n",
    "plt.xlabel('Time, s')\n",
    "plt.ylabel('Amplitude')\n",
    "plt.legend()\n",
    "plt.grid()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## We learned:\n",
    "\n",
    "* Approximation of numerical first derivative of a time dependent function by forward, backward and central FD operators\n",
    "* Comparison with analytical solution to estimate the numerical errors \n",
    "* Approximation of 2nd derivative of a time dependent function using Taylor series expansion"
   ]
  }
 ],
 "metadata": {
  "anaconda-cloud": {},
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}
